{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6I2vna9B7jXY"
      },
      "source": [
        "# Parcels trajectory code\n",
        "\n",
        "This notebook should provide a way to run the parcels code using local forcing fields.  The program will access data from the APDRC ERDDAP server, creat local netcdf files, read forcing from these files, then run trajectories in a forward mode.\n",
        "\n",
        "- JimP June 2022\n",
        "\n",
        "- Updated June 13 to run on google colab"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 1. Read/write to google drive\n",
        "\n",
        "To run on google colab we first need to mount the directory on google drive.  This will be used to read forcing files and also to write output and save the plots.  The block below is standard.  Later, to access, we use \"/content/drive/MyDrive/REU_2022\".  I've added sub folders for forcing, output and plots.\n"
      ],
      "metadata": {
        "id": "uCZpiWxXHlAe"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "KNJ05uoC8Q0N",
        "outputId": "6d5e837f-16c2-4500-f269-29f573373954"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Mounted at /content/drive\n"
          ]
        }
      ],
      "source": [
        "from google.colab import drive\n",
        "drive.mount('/content/drive')"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 2. Install \"missing\" packages\n",
        "\n",
        "Some additional packages are needed since they don't come with google colab; load these here using pip (need the leading !)"
      ],
      "metadata": {
        "id": "ajXT1giJIcAM"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 1000
        },
        "id": "ly1cL5vO7nTr",
        "outputId": "a3cb0fac-f6ec-4c26-c9a4-a87ab06901a2"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Looking in indexes: https://pypi.org/simple, https://us-python.pkg.dev/colab-wheels/public/simple/\n",
            "Collecting basemap\n",
            "  Downloading basemap-1.3.3-cp37-cp37m-manylinux1_x86_64.whl (863 kB)\n",
            "\u001b[K     |████████████████████████████████| 863 kB 4.5 MB/s \n",
            "\u001b[?25hRequirement already satisfied: numpy<1.23,>=1.21 in /usr/local/lib/python3.7/dist-packages (from basemap) (1.21.6)\n",
            "Requirement already satisfied: matplotlib<3.6,>=1.5 in /usr/local/lib/python3.7/dist-packages (from basemap) (3.2.2)\n",
            "Collecting pyproj<3.4.0,>=1.9.3\n",
            "  Downloading pyproj-3.2.1-cp37-cp37m-manylinux2010_x86_64.whl (6.3 MB)\n",
            "\u001b[K     |████████████████████████████████| 6.3 MB 32.7 MB/s \n",
            "\u001b[?25hCollecting basemap-data<1.4,>=1.3.2\n",
            "  Downloading basemap_data-1.3.2-py2.py3-none-any.whl (30.5 MB)\n",
            "\u001b[K     |████████████████████████████████| 30.5 MB 53.6 MB/s \n",
            "\u001b[?25hCollecting pyshp<2.2,>=1.2\n",
            "  Downloading pyshp-2.1.3.tar.gz (219 kB)\n",
            "\u001b[K     |████████████████████████████████| 219 kB 43.6 MB/s \n",
            "\u001b[?25hRequirement already satisfied: pyparsing!=2.0.4,!=2.1.2,!=2.1.6,>=2.0.1 in /usr/local/lib/python3.7/dist-packages (from matplotlib<3.6,>=1.5->basemap) (3.0.9)\n",
            "Requirement already satisfied: kiwisolver>=1.0.1 in /usr/local/lib/python3.7/dist-packages (from matplotlib<3.6,>=1.5->basemap) (1.4.4)\n",
            "Requirement already satisfied: python-dateutil>=2.1 in /usr/local/lib/python3.7/dist-packages (from matplotlib<3.6,>=1.5->basemap) (2.8.2)\n",
            "Requirement already satisfied: cycler>=0.10 in /usr/local/lib/python3.7/dist-packages (from matplotlib<3.6,>=1.5->basemap) (0.11.0)\n",
            "Requirement already satisfied: typing-extensions in /usr/local/lib/python3.7/dist-packages (from kiwisolver>=1.0.1->matplotlib<3.6,>=1.5->basemap) (4.1.1)\n",
            "Requirement already satisfied: certifi in /usr/local/lib/python3.7/dist-packages (from pyproj<3.4.0,>=1.9.3->basemap) (2022.6.15)\n",
            "Requirement already satisfied: six>=1.5 in /usr/local/lib/python3.7/dist-packages (from python-dateutil>=2.1->matplotlib<3.6,>=1.5->basemap) (1.15.0)\n",
            "Building wheels for collected packages: pyshp\n",
            "  Building wheel for pyshp (setup.py) ... \u001b[?25l\u001b[?25hdone\n",
            "  Created wheel for pyshp: filename=pyshp-2.1.3-py3-none-any.whl size=37324 sha256=edc33fcf521d1bed4599c72b5e4eed67b1ca444f5e2db6e04c940fa6667206b5\n",
            "  Stored in directory: /root/.cache/pip/wheels/43/f8/87/53c8cd41545ba20e536ea29a8fcb5431b5f477ca50d5dffbbe\n",
            "Successfully built pyshp\n",
            "Installing collected packages: pyshp, pyproj, basemap-data, basemap\n",
            "Successfully installed basemap-1.3.3 basemap-data-1.3.2 pyproj-3.2.1 pyshp-2.1.3\n"
          ]
        },
        {
          "output_type": "display_data",
          "data": {
            "application/vnd.colab-display-data+json": {
              "pip_warning": {
                "packages": [
                  "mpl_toolkits"
                ]
              }
            }
          },
          "metadata": {}
        },
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Looking in indexes: https://pypi.org/simple, https://us-python.pkg.dev/colab-wheels/public/simple/\n",
            "Collecting parcels\n",
            "  Downloading parcels-2.0.0.tar.gz (1.4 MB)\n",
            "\u001b[K     |████████████████████████████████| 1.4 MB 4.2 MB/s \n",
            "\u001b[?25hBuilding wheels for collected packages: parcels\n",
            "  Building wheel for parcels (setup.py) ... \u001b[?25l\u001b[?25hdone\n",
            "  Created wheel for parcels: filename=parcels-2.0.0-py3-none-any.whl size=1443240 sha256=65b919c19020f6232ecc2b361dc8b50ebcadb8ea7cde52a197474b40399527e3\n",
            "  Stored in directory: /root/.cache/pip/wheels/0a/24/2b/a319803019867c4186d9f91f4c71d094aefc5c1cec303f6cb9\n",
            "Successfully built parcels\n",
            "Installing collected packages: parcels\n",
            "Successfully installed parcels-2.0.0\n",
            "Looking in indexes: https://pypi.org/simple, https://us-python.pkg.dev/colab-wheels/public/simple/\n",
            "Collecting cgen\n",
            "  Downloading cgen-2020.1.tar.gz (19 kB)\n",
            "Collecting pytools>=2015.1.2\n",
            "  Downloading pytools-2022.1.12.tar.gz (70 kB)\n",
            "\u001b[K     |████████████████████████████████| 70 kB 3.1 MB/s \n",
            "\u001b[?25hRequirement already satisfied: numpy>=1.6 in /usr/local/lib/python3.7/dist-packages (from cgen) (1.21.6)\n",
            "Collecting platformdirs>=2.2.0\n",
            "  Downloading platformdirs-2.5.2-py3-none-any.whl (14 kB)\n",
            "Requirement already satisfied: typing_extensions>=4.0 in /usr/local/lib/python3.7/dist-packages (from pytools>=2015.1.2->cgen) (4.1.1)\n",
            "Building wheels for collected packages: cgen, pytools\n",
            "  Building wheel for cgen (setup.py) ... \u001b[?25l\u001b[?25hdone\n",
            "  Created wheel for cgen: filename=cgen-2020.1-py3-none-any.whl size=15836 sha256=aaa29a6ce27a707410693244d8d33cd67390073c8a3aa5304dab4acf1fa9186d\n",
            "  Stored in directory: /root/.cache/pip/wheels/2a/ea/03/babc69df7b9bbdd77352d860912ae9a246b538f3fef47792c5\n",
            "  Building wheel for pytools (setup.py) ... \u001b[?25l\u001b[?25hdone\n",
            "  Created wheel for pytools: filename=pytools-2022.1.12-py2.py3-none-any.whl size=65034 sha256=6fe2a2169ce7249c3a119fa7835020c470a6920a0bfdc7d6244f3e0723e98863\n",
            "  Stored in directory: /root/.cache/pip/wheels/37/5e/9e/76d7430e116b7cab0016fbabb26b896daae1946a3f7dea9915\n",
            "Successfully built cgen pytools\n",
            "Installing collected packages: platformdirs, pytools, cgen\n",
            "Successfully installed cgen-2020.1 platformdirs-2.5.2 pytools-2022.1.12\n",
            "Looking in indexes: https://pypi.org/simple, https://us-python.pkg.dev/colab-wheels/public/simple/\n",
            "Collecting pymbolic\n",
            "  Downloading pymbolic-2022.1.tar.gz (96 kB)\n",
            "\u001b[K     |████████████████████████████████| 96 kB 3.0 MB/s \n",
            "\u001b[?25hRequirement already satisfied: pytools>=2 in /usr/local/lib/python3.7/dist-packages (from pymbolic) (2022.1.12)\n",
            "Requirement already satisfied: pytest>=2.3 in /usr/local/lib/python3.7/dist-packages (from pymbolic) (3.6.4)\n",
            "Requirement already satisfied: six>=1.10.0 in /usr/local/lib/python3.7/dist-packages (from pytest>=2.3->pymbolic) (1.15.0)\n",
            "Requirement already satisfied: more-itertools>=4.0.0 in /usr/local/lib/python3.7/dist-packages (from pytest>=2.3->pymbolic) (8.13.0)\n",
            "Requirement already satisfied: pluggy<0.8,>=0.5 in /usr/local/lib/python3.7/dist-packages (from pytest>=2.3->pymbolic) (0.7.1)\n",
            "Requirement already satisfied: setuptools in /usr/local/lib/python3.7/dist-packages (from pytest>=2.3->pymbolic) (57.4.0)\n",
            "Requirement already satisfied: py>=1.5.0 in /usr/local/lib/python3.7/dist-packages (from pytest>=2.3->pymbolic) (1.11.0)\n",
            "Requirement already satisfied: attrs>=17.4.0 in /usr/local/lib/python3.7/dist-packages (from pytest>=2.3->pymbolic) (22.1.0)\n",
            "Requirement already satisfied: atomicwrites>=1.0 in /usr/local/lib/python3.7/dist-packages (from pytest>=2.3->pymbolic) (1.4.1)\n",
            "Requirement already satisfied: typing-extensions>=4.0 in /usr/local/lib/python3.7/dist-packages (from pytools>=2->pymbolic) (4.1.1)\n",
            "Requirement already satisfied: platformdirs>=2.2.0 in /usr/local/lib/python3.7/dist-packages (from pytools>=2->pymbolic) (2.5.2)\n",
            "Building wheels for collected packages: pymbolic\n",
            "  Building wheel for pymbolic (setup.py) ... \u001b[?25l\u001b[?25hdone\n",
            "  Created wheel for pymbolic: filename=pymbolic-2022.1-py3-none-any.whl size=118952 sha256=53cffa77013f15291520a5b5409f00d3e6f989e8c100563220f9c487b45def00\n",
            "  Stored in directory: /root/.cache/pip/wheels/e3/a0/06/2b011ece724068aa847ba6e0bfef76cec1f6060f098ffa81a9\n",
            "Successfully built pymbolic\n",
            "Installing collected packages: pymbolic\n",
            "Successfully installed pymbolic-2022.1\n"
          ]
        }
      ],
      "source": [
        "!pip install basemap\n",
        "!pip install parcels\n",
        "!pip install cgen\n",
        "!pip install pymbolic"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 3. Import all packages"
      ],
      "metadata": {
        "id": "VqCgi3E_I-mN"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "m-jhmYFk7jXq"
      },
      "outputs": [],
      "source": [
        "# load all necessary packages\n",
        "from parcels import FieldSet, ParticleSet, Variable, JITParticle, AdvectionRK4, plotTrajectoriesFile, ErrorCode\n",
        "import numpy as np\n",
        "import pandas as pd\n",
        "import math\n",
        "from datetime import timedelta\n",
        "from operator import attrgetter\n",
        "import pymbolic\n",
        "import matplotlib.pyplot as plt\n",
        "from mpl_toolkits.basemap import Basemap\n",
        "from netCDF4 import Dataset, num2date"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "0EEN7p0n7jXd"
      },
      "source": [
        "## 4. Set run parameters\n",
        "The variables \"model\" and \"depth\" will be used as the file name for the forcing and output.  A forward run will release particles along a line specified by lat1/lat2 and lon1/lon2 with \"npart\" number of points\n",
        "\n",
        "There are three basic steps that can be enabled:\n",
        "1. read forcing from server (read_data = True) or from local disk (read_data = False)\n",
        "2. run parcels model (run_parcels = True)\n",
        "3. make plots (make_plots = True)\n",
        "4. release multiple trajectories from same point (repeat = True)\n",
        "For example, to simply make plots from older runs, set run_parcels to False and make_plots to True"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "i9kMS8M-7jXc"
      },
      "source": [
        "### 4A. Note on forcing fields\n",
        "The forcing fields used come from two different global, operational models: the European Center for Medium Range Forecasts (ECMWF) Ocean Reanalysis (ORA-S5) and the NOAA Global Ocean Data Assimliation System (GODAS).  We get the output from the APDRC ERDDAP server as follows:\n",
        "* Monthly mean GODAS zonal (east-west) current: http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_fbd0_32c1_b107.html\n",
        "*  Monthly mean GODAS meridional (north-south) current: http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_3720_6ac2_d507.html\n",
        "* Monthly mean ORA-S5 zonal (east-west) current: http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_9faa_add2_cf6f.html\n",
        "* Monthly mean ORA-S5 meridional (north-south) current: http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_5e69_05b4_223f.html\n",
        "\n",
        "The process for creating and saving a local version:\n",
        "1. go to specific URL above\n",
        "2. select the time range\n",
        "3. select the lat range: 0 to 60N\n",
        "4. select the lon range: 120 to 280\n",
        "5. click output as netcdf (.nc), OR,\n",
        "5. click on \"generate the URL\"\n",
        "6. copy/paste the resulting URL in the notebook (below is an example)\n",
        "\n",
        "http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_fbd0_32c1_b107.htmlTable?uogrddsl[(2001-01-01T00:00:00Z):1:(2015-12-31)][(105):1:(105)][(0):1:(60)][(120):1:(280)]\n",
        "\n",
        "http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_3720_6ac2_d507.htmlTable?vogrddsl[(2001-01-01T00:00:00Z):1:(2015-12-31)][(105):1:(105)][(0):1:(60)][(120):1:(280)]\n",
        "\n",
        "http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_9faa_add2_cf6f.htmlTable?vozocrte[(2001-01-01):1:(2015-12-31)][(108.0303):1:(108.0303)][(0):1:(60)][(120):1:(280)]\n",
        "\n",
        "http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_5e69_05b4_223f.htmlTable?vomecrtn[(2001-01-01):1:(2015-12-31)][(108.0303):1:(108.0303)][(0):1:(60)][(120):1:(280)]\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "84fbqgnG7jXf"
      },
      "outputs": [],
      "source": [
        "read_data = False\n",
        "run_parcels = True\n",
        "make_plots = True\n",
        "repeat = False\n",
        "forward=True"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "usJuRRqO7jXg"
      },
      "source": [
        "The model start time is specified by year and month, and model can be ORA or GODAS"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "m5YCpUiJ7jXh"
      },
      "outputs": [],
      "source": [
        "year = '2022'\n",
        "month1 = '01'\n",
        "month2 = '12'\n",
        "model = 'ORA'\n",
        "vdepth = 5.140361 #1045.854 #508.6399 108.0303 1945.296 47.21189 1516.364 5.140361 26.5583 13.99104\n",
        "#model = 'GODAS'\n",
        "#vdepth =  949  #459.0 105.0 949 2174 55 1479 5 25 15\n",
        "depth = str(int(vdepth))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "CRPSYmqy7jXi"
      },
      "source": [
        "Set length of run in days"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "L8GvT67x7jXi"
      },
      "outputs": [],
      "source": [
        "yrsruntime= 25.0\n",
        "lengthofrun= str(int(yrsruntime)) # in years\n",
        "ndays = 365.0 * yrsruntime"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ncT5yhF_7jXj"
      },
      "source": [
        "Set number of repeat trajectories.  For example, if this is set to 24 hours, a trajectory will be released at this point once per day"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Jk9Ms8sn7jXk"
      },
      "outputs": [],
      "source": [
        "repeat_hour = 30.0 * 24.0\n",
        "nstart = 1\n",
        "if repeat == False:\n",
        "    repeat_hour = 0"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "uUUK35u47jXl"
      },
      "source": [
        "Setup parameters \n",
        "* specify the number of starting points (npart), starting lat and lon (lat1,lon1,lat2,lon2), OR\n",
        "* provide list of lats/lons for starting points"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "m7mgokix7jXm",
        "outputId": "4d8c91b3-898c-4dcd-a33a-703976c288db"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Number of latitudes 100\n",
            "Number of longitudes 100\n"
          ]
        }
      ],
      "source": [
        "if forward:\n",
        "  experiment = 'forward_' + model + '-' + year + '_' + depth + '_' + lengthofrun + '_years_' + 'box'\n",
        "else:\n",
        "  experiment = 'backward_' + model + '-' + year + '_' + depth + '_' + lengthofrun + '_years_' + 'box'\n",
        "\n",
        "'''\n",
        "# for line\n",
        "npart = 10\n",
        "lat1 = 30.0\n",
        "lat2 = 30.0\n",
        "lon1 = 210.0\n",
        "lon2 = 220.0\n",
        "'''\n",
        "# start lat and lon of garbage patch particles\n",
        "\n",
        "start_lat = [31.0, 31.0, 31.0, 31.0, 31.0, 31.0, 31.0, 31.0, 31.0, 31.0, \n",
        "             32.0, 32.0, 32.0, 32.0, 32.0, 32.0, 32.0, 32.0, 32.0, 32.0, \n",
        "             33.0, 33.0, 33.0, 33.0, 33.0, 33.0, 33.0, 33.0, 33.0, 33.0, \n",
        "             34.0, 34.0, 34.0, 34.0, 34.0, 34.0, 34.0, 34.0, 34.0, 34.0, \n",
        "             35.0, 35.0, 35.0, 35.0, 35.0, 35.0, 35.0, 35.0, 35.0, 35.0, \n",
        "             36.0, 36.0, 36.0, 36.0, 36.0, 36.0, 36.0, 36.0, 36.0, 36.0, \n",
        "             37.0, 37.0, 37.0, 37.0, 37.0, 37.0, 37.0, 37.0, 37.0, 37.0, \n",
        "             38.0, 38.0, 38.0, 38.0, 38.0, 38.0, 38.0, 38.0, 38.0, 38.0, \n",
        "             39.0, 39.0, 39.0, 39.0, 39.0, 39.0, 39.0, 39.0, 39.0, 39.0, \n",
        "             40.0, 40.0, 40.0, 40.0, 40.0, 40.0, 40.0, 40.0, 40.0, 40.0]\n",
        "\n",
        "start_lon = [211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220,\n",
        "             211, 212, 213,214,215,216,217,218,219,220]\n",
        "\n",
        "'''\n",
        "# start lat and lon of Mexico particles\n",
        "\n",
        "start_lat=[22,\n",
        "           21,\n",
        "           20,\n",
        "           19,\n",
        "           17,\n",
        "           15,\n",
        "           15,\n",
        "           14,\n",
        "           15,\n",
        "           13,\n",
        "           12,\n",
        "           10,\n",
        "           8,\n",
        "           7,\n",
        "           7,\n",
        "           7]\n",
        "start_lon=[250,252,254,256,258,260,262,264,266,268,270,272,274,276,278,280]\n",
        "\n",
        "\n",
        "# start lat and lon of Japan particles\n",
        "start_lat=[30,\n",
        "           31,\n",
        "           31,\n",
        "           32,\n",
        "           33,\n",
        "           32, \n",
        "           32,\n",
        "           32,\n",
        "           32,\n",
        "           32,\n",
        "           33,\n",
        "           34,\n",
        "           36,\n",
        "           38, \n",
        "           40,\n",
        "           40]\n",
        "\n",
        "start_lon=[130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 140.2, 140,4, 140.6, 143.8]\n",
        "\n",
        "\n",
        "# start lat and lon of Philippines particles\n",
        "start_lat=\n",
        "'''\n",
        "\n",
        "print('Number of latitudes', len(start_lat))\n",
        "print('Number of longitudes',len(start_lon))\n",
        "\n",
        "npart=int(len(start_lat))\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "fDP2ks3J7jXn",
        "outputId": "143e94b7-7adf-4b8f-95f2-dfb1ce7a524f"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Running experiment  forward_ORA-2022_5_25_years_box\n"
          ]
        }
      ],
      "source": [
        "print('Running experiment ', experiment)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Mv_Wwh3G7jXn"
      },
      "source": [
        "### This should be the end of the user-defined parameters...."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "IBCxXAta7jXo"
      },
      "source": [
        "## 5. Get forcing files\n",
        "\n",
        "If the variable \"read_data\" is true, then we go out to the server to get the data at run time.  Otherwise, this step is skipped and we read directly from disk (google drive)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "EXcBxUwN7jXo"
      },
      "outputs": [],
      "source": [
        "if read_data:\n",
        "    import urllib\n",
        "    time1 = year + '-' + month + '-' + day1\n",
        "    time2 = year + '-' + month + '-' + day2\n",
        "    constraint = '[(' + time1 + '):1:(' + time2 + ')][(vdepth):1:(vdepth)][(0):1:(60)][(120):1:(280)]'\n",
        "    if model == 'ORA':\n",
        "        root_url = 'http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_'\n",
        "        URL_U = root_url + '9faa_add2_cf6f.htmlTable?vozocrte' + constraint\n",
        "        URL_V = root_url + '5e69_05b4_223f.htmlTable?vomecrtn' + constraint\n",
        "    elif model == 'GODAS':\n",
        "        root_url = 'http://apdrc.soest.hawaii.edu/erddap/griddap/hawaii_soest_'\n",
        "        URL_U = root_url + 'fbd0_32c1_b107.htmlTable?uogrddsl' + constraint\n",
        "        URL_V = root_url + '3720_6ac2_d507.htmlTable?vogrddsl' + constraint\n",
        "    \n",
        "    urllib.request.urlretrieve(URL_U, '/content/drive/MyDrive/REU_2022_copy/forcing_copy/' + model + '-' + year + '_' + depth + '_U_forcing.nc')\n",
        "    urllib.request.urlretrieve(URL_V, '/content/drive/MyDrive/REU_2022_copy/forcing_copy/' + model + '-' + year + '_' + depth + '_V_forcing.nc')"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "l_fKrn577jXp"
      },
      "outputs": [],
      "source": [
        "# specify the filenames for the U,V fields\n",
        "filenames = {'U': '/content/drive/MyDrive/REU_2022_copy/forcing_copy/' + model + '-' + year + '_' + depth + '_U_forcing.nc',\n",
        "             'V': '/content/drive/MyDrive/REU_2022_copy/forcing_copy/' + model + '-' + year + '_' + depth + '_V_forcing.nc'}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "YQaEDS437jXp"
      },
      "outputs": [],
      "source": [
        "# set the variable names; here the first part is what parcels is expecting, the second is from the file\n",
        "#   'U' is the zonal velocity, given in our file by 'vozocrte'\n",
        "#   'V' is the meridional velocity, given in our file by 'vomecrtn'\n",
        "#   'lat' is the latitude, given in our file by 'latitude'\n",
        "#   'lon' is the longitude, given in our file by 'longitude'\n",
        "#   'time' is time, given in our file by 'time'\n",
        "if model == 'ORA':\n",
        "    variables = {'U': 'vozocrte',\n",
        "                 'V': 'vomecrtn'}\n",
        "if model == 'GODAS':\n",
        "    variables = {'U': 'uogrddsl',\n",
        "                 'V': 'vogrddsl'}\n",
        "\n",
        "dimensions = {'lat': 'latitude',\n",
        "              'lon': 'longitude',\n",
        "              'time': 'time'}"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "y_ODX7Su7jXq"
      },
      "source": [
        "## 6. Initialize parcels run"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "FWxTJ7UB7jXq",
        "outputId": "ce126c18-624d-41cd-f4a4-153770990a52"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Forcing fields (u,v) start at 1979-01-15 00:00:00 and 1979-01-15 00:00:00\n",
            "Forcing fields (u,v) end at 2018-12-15 00:00:00 and 2018-12-15 00:00:00\n"
          ]
        }
      ],
      "source": [
        "# Check start and end times\n",
        "uvar = Dataset(filenames['U'])\n",
        "vvar = Dataset(filenames['V'])\n",
        "timevar = uvar.variables['time']\n",
        "time_convert = np.array(num2date(timevar[:],timevar.units))\n",
        "time_convert = time_convert.astype('datetime64[ns]')\n",
        "udate1 = pd.to_datetime(time_convert)[0]\n",
        "udate2 = pd.to_datetime(time_convert)[-1]\n",
        "timevar = vvar.variables['time']\n",
        "time_convert = np.array(num2date(timevar[:],timevar.units))\n",
        "time_convert = time_convert.astype('datetime64[ns]')\n",
        "vdate1 = pd.to_datetime(time_convert)[0]\n",
        "vdate2 = pd.to_datetime(time_convert)[-1]\n",
        "print('Forcing fields (u,v) start at', udate1, 'and', vdate1)\n",
        "print('Forcing fields (u,v) end at', udate2, 'and', vdate2)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Strt6a197jXr"
      },
      "outputs": [],
      "source": [
        "# create fieldset and run forward\n",
        "if run_parcels:\n",
        "    fieldset = FieldSet.from_netcdf(filenames, variables, dimensions,allow_time_extrapolation=True) # check to see if allow_time_extrapolation only refers to the last month or does it go through all 12 months after\n",
        "    # time-step in forcing file in seconds\n",
        "    dt_forcing = fieldset.U.grid.time[1] - fieldset.U.grid.time[0] # 2678400 seconds = 31 days = 1 month\n",
        "    repeatdt = timedelta(hours=repeat_hour)\n",
        "    \n",
        "    # start_lon = np.linspace(lon1,lon2,npart) # use only if start_lat and start_lon are not giving specific points\n",
        "    # start_lat = np.linspace(lat1,lat2,npart)\n",
        "    if repeat:\n",
        "       #from a line\n",
        "       # pset = ParticleSet.from_line(fieldset=fieldset, pclass=JITParticle, size=npart,\n",
        "       #                              start=(lon1,lat1), finish=(lon2,lat2), #time= yrsruntime*12*2678400, # runtime (years) in seconds\n",
        "       #                              repeatdt=repeatdt)\n",
        "        #group of points\n",
        "        pset = ParticleSet(fieldset=fieldset, pclass=JITParticle, \n",
        "                                     lon=start_lon, lat=start_lat, depth=None, #time= yrsruntime*12*2678400,\n",
        "                                     repeatdt=repeatdt)\n",
        "\n",
        "                                   \n",
        "    else:\n",
        "      # from a line\n",
        "      #pset = ParticleSet.from_line(fieldset=fieldset, pclass=JITParticle, size=npart,\n",
        "      #                            start=(lon1,lat1), finish=(lon2,lat2), # time= yrsruntime*12*2678400,\n",
        "      #                            repeatdt=None)\n",
        "      # group of points\n",
        "      pset = ParticleSet(fieldset=fieldset, pclass=JITParticle, \n",
        "                                     lon=start_lon, lat=start_lat, depth=None,\n",
        "                                     repeatdt=None)\n",
        "  \n",
        "\n",
        "                                    \n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "gJha4X4H7jXr",
        "outputId": "c8c1eec8-d2b9-4974-d30a-8c37261df6ca"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Run should have trajectories starting/ending at 100 locations with 1 start times at each and 295 integration (path) times\n",
            "Arrays will be 100 by 295\n"
          ]
        }
      ],
      "source": [
        "months= 12*(yrsruntime)\n",
        "ntime = months / ( dt_forcing / 86400.0 * months ) * ndays + 1\n",
        "#ntime = 24.0 / ( dt_forcing / 86400.0 * 24 ) * ndays + 1 #what is the purpose of this?\n",
        "ntraj = npart\n",
        "if repeat:\n",
        "    nstart = ( ndays * 24 / repeat_hour + 1 )\n",
        "print('Run should have trajectories starting/ending at', int(ntraj), 'locations with', int(nstart), \n",
        "      'start times at each and', int(ntime), 'integration (path) times')\n",
        "print('Arrays will be',int(ntraj*nstart),'by',int(ntime))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "fnEAVK127jXs"
      },
      "source": [
        "## 7. Run code\n",
        "__NOTE__: sometime you have to run this twice (?) on local machines; so far this is NOT The case on google colab\n",
        "\n",
        "The output file contains the lat/lon of the particles over time (set here as traj_output)."
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "This is taken from [tutorial_Agulhasparticles](https://github.com/OceanParcels/parcels/blob/master/parcels/examples/tutorial_Agulhasparticles.ipynb) in order to delete particles that go out of bounds. I don't know if there is another option that can be used to track the particles after they leave the perimeter because there is no data in the netCDF file that has velocities for those positions. "
      ],
      "metadata": {
        "id": "8MFuWbZZsHVp"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "def DeleteParticle(particle,fieldset,time):\n",
        "  particle.delete()"
      ],
      "metadata": {
        "id": "J4Dg6QE4siMm"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "MZaFL5Np7jXs",
        "outputId": "ac1a7fac-fbd2-4db2-85e4-3ad4fa65512a"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "INFO: Compiled JITParticleAdvectionRK4 ==> /tmp/parcels-0/6a1e0248a8d3cd95183d439bd993946e.so\n",
            "100% (788400000.0 of 788400000.0) |######| Elapsed Time: 0:00:40 Time:  0:00:40\n"
          ]
        }
      ],
      "source": [
        "# specify the output file that will contain the lat/lon points of the particles then integrate\n",
        "# here we set the output time resolution to be the same as the forcing (e.g., tides at 1 hour\n",
        "# and ROMS at 3 hours)\n",
        "\n",
        "traj_output = '/content/drive/MyDrive/REU_2022_copy/output_copy/box_output_copy/' + experiment + '.nc'\n",
        "if run_parcels:\n",
        "  if forward:\n",
        "    output_file = pset.ParticleFile(name=traj_output, outputdt=timedelta(hours=dt_forcing/3600.0)) #dt_forcing/3600=output every month\n",
        "    pset.execute(AdvectionRK4, runtime=timedelta(days=ndays), \n",
        "                 dt=timedelta(hours=24), \n",
        "                 recovery={ErrorCode.ErrorOutOfBounds:DeleteParticle},\n",
        "                 output_file=output_file)\n",
        "# run backwards by making timedelta -timedelta and then make time=the last time in forcing field when we define the pset\n",
        "  else:\n",
        "    output_file = pset.ParticleFile(name=traj_output, outputdt=timedelta(hours=dt_forcing/3600.0)) #dt_forcing/3600=output every month\n",
        "    pset.execute(AdvectionRK4, runtime=timedelta(days=ndays), \n",
        "                 dt=-(timedelta(hours=24)), \n",
        "                 recovery={ErrorCode.ErrorOutOfBounds:DeleteParticle},\n",
        "                 output_file=output_file)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "x-0nMgy07jXs"
      },
      "source": [
        "## 8. Plot results"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "U7b881sW7jXt",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 326
        },
        "outputId": "83a6d802-4e1d-4e11-cd59-06f20226b8a9"
      },
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 720x720 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlsAAAE1CAYAAADQ0yXGAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOydd3gUVfuw75MKCRAIkEBC74TQQQEBgQCK0qSFDlEUpQkoij8LChZeEQJopChKlSI1IL2FKr1J74SeBBJCEtL2fH+kfAkzKZu2u+Hc18W1ZGbOzLOwmX3OPc85R0gpUSgUCoVCoVDkDlamDkChUCgUCoUiP6OSLYVCoVAoFIpcRCVbCoVCoVAoFLmISrYUCoVCoVAochGVbCkUCoVCoVDkIirZUigUCoVCochFVLKlUCgUCoVCkYuoZEuhUGQKIcRgIcQZIUSkEOK+EGKWEKJoiv1fCyFihRBPhRChQogDQoimOueZL4SIE0KUzuB6Kc+X9KdSDr2XQUKIY0KIJ0KI20KIH4UQNin27xZCPEtx3Ys5cV2FQvFiopIthUKRIUKIj4D/AeMAJ6AJUB7YJoSwS3HocillIaAEsAv4+7nzOALdgTCgfyYuvVxKWSjFn2vZfzcAOACjE+N8GfACPn7umBEprls9h66bZVImgwqFwrJQyZZCoUgXIUQR4BtgpJRys5QyVkp5A+gFVEAnaZJSxgFLAHchRMkUu7oDocBEYFA2YqoghJBCCB8hRKAQ4rEQ4n0hRGMhxOlEs/ZLWu2llLOklHullDFSyjuJsb6ShThKJZq+4im2NRBCBAkhbBN/flsIcT4xxi1CiPIpjp2RGP+TRNPWIsW+r4UQK4UQi4UQT4DBQoiXhBBHE49/IISYZmzMCoUi71HJlkKhyIhmQAFgdcqNUsqnwEag3fMNEm3XQCAEeJxi1yBgKbAMqCGEaJjBtTsJIR4JIc4KIT7Q2f8yUBXwBqYDnwNtgVpALyHEqxm/PQBaAmef2/aDECJYCLFfCNFKr5GU8j6wm4TEM4kBwDIpZawQogvwf0A3oCSwl4T3n8QRoB7gDPwF/C2EKJBifxdgJVCUhIRwBjBDSlkEqAysyOT7UygUJkQlWwqFIiNKAMGJtup57iXuT6KXECIUiALeBXoktRNClANaA39JKR8AO0hIyNJiBVCThCTlXeArIUSf546ZJKV8JqXcCkQAS6WUDxNt1V6gfkZvTgjxNtAI+CnF5k+BSoA7MBdYL4SonMYpFpBo94QQ1kAfYFHivveBH6SU5xP/Hb4H6iXZLSnlYilliJQyTko5FbAHUj6yPCilXCulNEgpo4BYoIoQooSU8qmU8t+M3p9CoTA9KtlSKBQZEQyUSKNmqHTi/iRWSCmLAq7Af0BKczUAOC+lPJn48xKgb9LjtueRUp6TUt6VUsZLKQ+QYHV6PHfYgxR/j9L5uVB6b0wI0RX4AeggpUx+H1LKQ1LKcClltJRyAbAfeCON06wDPIQQFUmwfGFSysOJ+8oDMxIfa4YCjwBBQhKHEOLjxEeMYYn7nUidvAY+d613gGrABSHEESFEx/Ten0KhMA9UwaVCociIg0A0CY/Ckh9bCSEKAR1IeEyWCillsBDiPeCoEOIvKeU9EixWOSHE/cTDbIDiJCQx6zIRhyQhUckRhBCvA78Bb0opz2T12lLKZ0KIFSTYrRr8f6sFCcnSd1LKJTrXbwF8QkJx/lkppUEI8fi568jnrnUZ6COEsCLh/2OlEKK4lDIig/gVCoUJUWZLoVCki5QyjIQC+Z+FEK8LIWyFEBVISLxukzq5SNnuIrAF+CRxCojKwEsk1CjVAzxJqFPSfZQohOgihCgmEngJGEXmkrIMEUK0IcGsdU9hoZL2FRVCvCaEKCCEsBFC9COhpmtzOqdcCAwGOpP632M28JkQolbiuZ2EED0T9xUG4oAgwEYI8RVQJIO4+wshSkopDSQMNAAwZPyOFQqFKVFmS6FQZIiU8kchRAgJdU2VgSfAWqCflDI6naZTgJ0kPB5b97xBEkLMAPYKIZyllI+ea9sb+IOEOqbbwP8SH+nlBF8mxrRRiGSRtFdK2QGwBb4lwVLFAxeArlLKS2mdTEq5XwhhAI5LKW+m2L4m0QAuS6zTCgO2kTAlxhYSErhLJNSb+aJ9bPg8rwPThBAOwE2gd2Itl0KhMGOElDLjoxQKhUKRLkKInSQU//9u6lgUCoV5oZIthUKhyCZCiMYkGKuyUspwU8ejUCjMC1WzpVAoFNlACLEA2A6MVomWQqHQQ5kthUKhUCgUilxEmS2FQqFQKBSKXEQlWwqFQqFQKBS5iNlO/SCEUM83FQqFQqFQWAxSSt3Jj8022QLYvXu3qUMwK27evImzszOFCxc2WQzHjx/H3t4eDw8PUsxPZBT+/v5MnTo1x2I6dOgQjRs3xsrKivXr13PixAlefTXt9YcXLVrE3LlzsbJKELuRkZGMHDmSqVOnUrRo0eTj/v77b3r27JnWacya8PBwbt68iaenp6lDyRX27t2Lk5MTtWvX1nwOHz58yOTJk+nSpYumXVBQEAaDAVdX1xyP6ejRo9jb2zNixIhU27t160apUqWwsrIiNDSUd999N3nfxo0bWbVqFVeuXMnxeDIiNjaW7777jipVqvDPP//QokUL1q1bx6+//sro0aP5+OOP8zwmc+PcuXPUqFEj+V6hSLhf3r9/n0qVKmXrPBEREZw6dYrPPvsshyIzPel9J6pPkAXh7+/PrVu3THb9yMhIbty4gaurKytXrmTVqlXcuXPHqHM8e/Ysx5PFcePGERsby82bN1m/fn26iVZgYCB169ZNdfN0cHBg3rx5qRItwGITLUhIzH19fU0dRq7RokUL6tSpo3tzc3FxITpaf57V48ePc+DAgRyP58KFC4SHh2sSLYDVq1fz66+/8tprr1GrVq3k7UFBQVy8eJGVK1fmeDyZwdbWlq+//pqSJUtSvHhxtm/fTs+ePRk4cCC9e/c2SUzmxuzZs4mL01t//cXl/v37OfKZffLkCSVKlMj4wHyC2Y5GFEJIZbZSk5Nm68mTJ4SFhVG2bFkAdu3aReHChalUqRLOzs66bdavX8+wYcOSezSPHz9myZIlnD17FicnJ1555RWKFEl3tRECAwNZtWoV1apVIyoqigEDBtC4ceNsvZdDhw5Rt25dhg4dSr9+/bC3t0/z2KVLlzJ16lQcHR2zdU1zJ7+brYyYNm0aZcuWxcXFJdX23DBbt27d4vjx40ybNi3Nnm14eDijRo1i8ODBBAYGsmfPHsqVK8fw4cMpVqxYjsWSVc6fP88333xD165diY2NpVy5cka1f/jwITExMZQpUyaXIjQNymxpySmzdeHCBVxdXenatWsORWZ6hBCW+RhRkRp/f3/atGmTqnecFgaDgaCgIO7evUtQUBCPHj3C2toaa2trbGxsKFmyJKGhoTx58oSgoCDq1KmDi4sLhw8fZu/evcTGxhIfH0+pUqWoVKkSJUuWJCYmJtUvWLFixZJ78leuXGHx4sXcuXOHSpUq0bhxY2xtbTVxubu7884771CoUCHi4+NZtmxZtpOtcePG0bhxYzp06JBuopXUk8rviRb8f7M1b948U4diEry9vZk+fTodO3ZMtf348eNERkby1ltv5di1du/ezZw5c9J9hPDNN99Qs2ZNFi1aRO3atZk+fToODg45FkN2mTNnDh988AFr1qyhcuXKRidbkZGRHDx4kD59+uRShKZh9uzZ/PTTT9jZ2Zk6FLMhyWx98skn2TpPeHj4C9UZVMmWBdG5c+c0rRPA06dPWb16NU5OTtjZ2VGuXDk8PT2pWrUqZcqUwdraOtXxUkoGDhxI+fLl6datGwDNmzdP3h8XF8eFCxc4evQoW7duZdSoUWleu1SpUrRt25ZDhw6xfv16Dh06xJgxYzTHWVlZUahQIR49esT169c5deqUsf8MGtq1a4etrS2lS5dO97idO3e+MHUo5cuX1/33f1Fwd3fn6dOnmu0NGjTAYEi9bnNcXBwRERE4OTll6VpOTk7Y2KR9Kz148CABAQEMGzaMMWPG6HZCTM3QoUMpXbo0LVu2xM/Pj02bNvH6669nui6zRIkSHDt2LN8lW++//366/7cvIqVKlaJHjx7ZPk9ERMQL9RhRfYosiPTMVnx8PMuXL8fX15fixYtn6nxCCKZPn57moz8bGxs8PT0z1fsYOXIk1apVo3r16nzzzTfJ20NDQ7l+/Tp3797lyZMn2NnZYWtrS9myZWnYsCGDBw/OVKxpcf78eWbNmsWCBemvT5xk6tzd3bN1PUvhRTdbBoOByMjI5J/Pnj3LkSNHCA8Pp0yZMqnM1qlTpzhx4gRDhgxJdY7jx4/j7u6e4SPHkJAQnjx5kubvUYMGDTh06JBZP4qaM2cO3t7eNG3alBEjRrBr1y4WL16Mt7d3KqsTGRnJ2rVr6dWrF8HBwdy7d4/69evj6OjI1atXiYyMNCtjl12U2dKSU2YrLi6O2NjYHIrK/FE1WxZEejVbK1eu5L333qN27domiAyOHDnC77//jrOzM48fP05OqkqVKkW9evWSH1NmdQSjHpGRkQwdOpSXX36ZWrVqpftltmPHDrp27UqDBg1y7PrmzItes+Xr60vhwoUJDw/nv//+o23btvTq1Yuff/6ZS5cu4e3tnXzsihUriIuLw9vbG2tra4KCgti4cSMnTpxgypQpGZqoGzducOTIEQDatGmTbIktifPnz1O6dOlUg0Ru3rzJl19+SZcuXVIZiJUrV3Ls2DG8vLwIDAzEx8cHgC+++IJOnTrx8ssv53n8uYWq2dKSUzVb0dHRrFu3Dj8/vxyKzPSomq18Qlpma9++fbRs2dJkiRZA48aNkwtr3d3dczSp0kNKyfjx4+natSsTJ05Mt/cppeT+/fsvTKIFL7bZOn/+PMuWLaN27dp06dKF0aNHJ38enZ2dkVKyatUqunXrhhACIQSdO3fm1KlT3L17F3t7e958802qVauWqUd+FSpUoEKFCty4cYPz58/n9tvLFVKarSTKly/P7NmzGTduHJ6ennh4eAAJU1k8ePCA//77L9Xvee3atbl27Vq+SraU2dKSU2bL3t4eFxcX9u7dy/LlyxFCYGNjg5WVFQMHDqRu3bo5FLF5oMyWBaFntq5cucLdu3f58ssvTRhZ3jN37lyioqKoV69ehr3PkydPUqFCBd15l/IrL7LZ+uuvv3Bzc6NVq1aafXfu3CE+Pp47d+4wb948evXqxc6dO5k0aRI+Pj6MHz+emjVrMmTIEAYNGpTpTkNMTAyLFy9m3rx5Flnjo2e2kpBSMm3aNEJDQ2nbti2Q8Plau3YtTk5OdOrUCSEEK1aswMrKKlU9T0REhEUPSFFmS0tOmS1IeNy/cOFCHB0dk6faiYuLY/PmzUyZMiXb589r0jNb6hNkQTw/z1ZoaCiHDh3i888/N2FUec+xY8c4f/489erVAzKeC+fMmTN06tQpr8IzC/L7PFvp0bdvX91ECxIeJ69fv56mTZsyduxYfvrpJ5o2bYqtrS2LFy/G09OT3377jebNmxtlZ9euXctXX31lkYkWJJittKycEIKPPvqIpk2bsmTJEuLi4ihcuDBNmjRh69atyfckg8FAuXLlePDgAZDwpWzpn0E1z5aWnJpnCxIGTDVu3JiHDx8mb7OxsdEd3GLpKLNlQaQ0W3FxcSxYsAA/Pz+Tziif14SFhTFy5Eh8fHySe5vp9T5v3bpFWFiY7mST+ZkX2WylR5LZSpra4MGDBxQtWjTVlCEff/yxUcn56dOnKVCgQKqZ4S2N9MxWSq5cucI333zDyy+/zIkTJ/D09CQkJISQkBAuXbqEn58fc+fOxd7enqtXr1K+fHnN9BuWhDJbWnLSbKU8Z8qBFVu2bOGdd96hcuXKOXaNvECZrXxCSrO1YsUKvvrqqxcq0ZJSMm7cOHr06JHq5pde73Pfvn28/fbbeRWi2fAim630SDJbSbi6umrmZrOxsdFMD5EWZ8+e5ezZs5qRjJZGemYrJVWqVOHXX3+lSJEi/Prrr4wZM4Zvv/2WWbNm4efnR6VKlZg2bRq9evXi8ePHtGzZMg+izz2U2dKSk2YriedHsNavXz/V72l+wDKd9wtK0jxb27dvp2vXrlSpUsXUIeUJ0dHRHD16lKVLl9KwYUPNfEhpzYUTExPDvXv3ePbsWb4ajp4ZXvR5ttLCy8uL+Pj4dI8pXrw4T58+TXMqB4PBwOHDh7l8+TJt27bFz88v1weE5DZJ82xlhsKFCzNw4EDN9ho1aiT/vX79+gwYMAB/f39q1aqlWSLLUlDzbGnJqXm20sPFxYU9e/bk6jXyGsv79L/A+Pv7s3v3booVK0aHDh1MHU6uEhAQwPjx4xk5ciRjx45l//79vP7666lu6Emk1fu0s7OjZ8+eTJw4keHDh+Pn50dwcHBehG9ylNnS53mzpYeLiwthYWGa7VJKtm/fzrJly2jUqBG//fYbvXv31kwWbIlk1mwZw8CBA5k9ezaXL1/mu+++y9Fz5xXKbGnJDbOlh8Fg4NmzZ7l+nbxC1WxZEMePH+fChQvMnj3b4nvS6REVFcWwYcPo379/pnqVma2ruHPnDr///jvTp0+nTp06ORWuWaJqtvR5vmZLj82bN3P16tVU/3ZSSlasWIG3t3e6C51bKpmt2TKWkydPUr9+faZPn548oMWSUDVbWnKjZkuPU6dOUbVqVYsSC6pmKx8QHR3N3Llz6dOnT75OtAAWL15My5YtM63vM9v7DA0N5bXXXsv3iRYos5UWmTFbzs7OPHnyJPlnKSUrV66kd+/e+TLRgtwxWwD16tVj5MiRjB49mgkTJuT4+XMbZba05JXZqlWrFjt37sz16+QV6mG0BSClZNmyZfj5+Vlk79BYDh8+TP/+/TN9fGbqKm7fvs21a9cscu6WrKBqtvTJqGbr6NGj/PLLL8lz/gCsXr2aHj16WHyxd3oYU7NlLDNmzCAkJMSiDEUSqmZLS17UbEHCQJXw8HDN9tjYWO7fv0/ZsmVzPYacRJktC2DTpk0MHjwYf39/i52hOrMcOnQo3Uc8emSm9xkWFkbTpk3zvRVMQpktfdIyW1JKfvnlF5YuXYqPj09ycfyaNWvo3LlzmvN25Rdyy2wZDAbGjx+fvLKEpaHMlpa8MluxsbGpEt0zZ84wfvx4RowYwZQpU/jmm2+IiYnJ9ThyClWzZeacOHGCggULMnTo0FyrqzAnPvzwQzp16pSpZVKSyExdhcFgYM2aNcycOTMnwjR7VM2WPno1W6GhoYwfP546deqkWgpr3bp1tG/fntdff90UoeYpuXVvkVIyceJEbG1teeWVV3L03HmBqtnSklc1W9u2baN169acPn2aCxcu4OrqSrNmzShYsCCQ8LRiy5YtjBkzhvr16+dqLJlF1WxZKHfu3OHOnTsMHToUyL3ep7kQHByMlNKoRAsy1/u0srIiJibGonpC2UGZLX2eN1v//vsvo0eP5s0330yVaG3YsAEvL68XItGC3Lu3CCGYMGEC5cqVY+3atZhr5z4tlNnSkldm67///mPdunW4ubnRt29fvLy8khMtgDJlyjB48GD++usvvvvuO2JjY3M9puygzJaZEhkZyfLly/ntt9+Sk4/8brYmT55MlSpVKFmypFHtMtv7PH36NBUrVrToGa0zizJb+iSZrbJlyzJ9+nTu379Phw4dUj1e3rhxI02bNqVr164mjDRvyYt7y6FDh5g1axZ9+/bVTCRrriizpSWvzJYx3Lp1i+3bt/PRRx+ZdACUMlsWhsFgYNmyZfz444+pLE9+N1t37twxOtGCzPc+a9Wqxfbt27MSmsWhzJY+O3bsYNmyZQwdOhRHR0feeOONVInW5s2beemll16oRAvy5t7y8ssv891337Fo0SKCgoJy9Vo5hTJbWvLKbBlDuXLlGDRoEAsWLOCHH34wy/8zlWyZGVJKVq1axciRI3FxcUm1b+jQodSsWdNEkeUuUsos/4JkdsSQtbU1UVFRZq+bcwI1GlEfR0dHDh8+TOfOnalevXqqfdu2baN+/fp0797dRNGZjry6t7i7uzN37lx27drFxYsXc/162UWNRtSSV6MRjcXa2prOnTvj7u7OO++8w9mzZ00dUipUsmVGxMfHs3TpUrp3706jRo00+/Oz2QoKCqJQoUJZamtM77NmzZr5au6WtFBmS8v27dv57bffcHd316wpun37djw8POjVq5eJojMteXlvKViwIL/88gshISH8+++/eXLNrKLMlhZzNFspKV++PAMHDuT333/np59+ynB5rrxC1WyZCdHR0SxZsoSPP/6Y2rVr6x6Tn2u2du/ezfHjx2nYsKHRbdOrq7h58yYHDx7EysoKGxsbrK2t6dKlS76dnDIJVbOl5e+//+bGjRuUL18eV1fX5O27du2iUqVKuuv9vSiY4t5y7NgxVq5cSa1atRBCUKZMmTy7dmZRNVtazLFmKy2uXbvGnj17GD9+vMZk5waqZsvMCQ8PZ+HChXz77bdpJlqQv83WmTNnjJ5fK4n0ep+nT5/m66+/xs/PjxkzZjBt2rR8n2iBMlshISGEhISk2laqVCmOHz/OgQMHkrcFBAQk94RfZPL63mIwGPj5558pUqQI586dY8eOHcn7jh07pjuZpSlQZkuLuZutlFSqVIn+/fsza9YsfH19MRgMJotFPYw2MQ8fPmTDhg388ssvGfYqc3OWZ1Nz69atdBPN9EivrsLGxsakv2Cm4kWu2QoODmbUqFG4ubkRERGBra0ttra2PH36lMePH3P79m3atm3LyZMncXNzw8fHx9Qhm5y8vrf8/vvvBAYG0q5dO/r168eKFSs4ffo0derU4bPPPmP16tV5Fkt6qJotLeZas5UWNjY2vPXWW1y5coV33nmHzz//nCpVquR5HMpsmZDr16+zc+dO5syZkyl9n5/NVkxMTJZnd0+v92lra/vCzK2VkhfVbMXFxTFu3Dj69etH2bJl6d27N927d6dly5ZERETg5uZG27ZtmTFjBsWLF2fIkCGmDtksyMt7S1hYGL///jvTpk2jX79+ADRr1oy7d+8CCfeCrNZv5jTKbGmxJLOVkipVqtCvXz9mzJjBzz//nOedcJWym4jTp09z+/Zt/Pz8Ml0PkJ/NVnZGCKbX+7S2tiY6OjrL57ZUXlSz9fnnn+Pl5UVMTAxbt24lNDSUuLg4IiMj+eWXX4iIiCA+Pp5OnTrRoEEDU4drNuTlvaVAgQLs2LEj1SCF+/fv4+TklCfXNwZltrRYmtlKia2tLT169ODSpUu88847fPXVV1SsWDFPrq3MlgnYt28f0dHRTJ482ajCy/xqtp4+fZqtG1p6vU8bGxtltl4QRowYQWhoKGXLliUgIIA//viDEiVK0L59e6ZOnUrRokWTZ5BXiVZq8vLeYm9vrxkNeufOHYoVK8atW7fMagCQMltaLNVspaRatWr07duXqVOn8uuvvxIcHJzrnXKVsucxmzZtombNmgwaNMjotvnVbAUGBlK8ePEst8/IbL2IydaLZrb27dvHgQMHaNeuHQsWLMDBwQFXV1eGDx+e6jgvLy+zGQpuTpj63nLnzh2cnJxYtGgRb7/9tsnieB5ltrRYstlKiZ2dHT179uTy5ctMnjyZqKgo4uPjsbKySv4jhEj+ExcXx9dff53l7ypltvIIKSUrVqzg1VdfzVKiBfnXbJUpU4ZHjx5luX16vc8XNdl6kczWxYsXadGiBWXLlqVRo0bMnTuXyZMn6x77/NqIigRMfW+5d+8eO3bsYOrUqdjZ2TFjxgyzMErKbGnJD2YrJVWrVqVTp0706tWLPn364O3tTc+ePenevTvdunXjrbfeomvXrri4uPDgwYMsX0el7HlEfHw8jx8/pk2bNlk+h6l7n7lF4cKFc61my8bGRtVs5XOaNWvG9evXKVeuHPHx8djY2KRZYK3Mlj6mvrc8fPiQjh07UrZsWZ4+fUrLli3Nwigps6Ulv5itvEaZrTzCxsaGrl27MmbMGG7dusXu3bv57bffmDhxYqbPYereZ26Scg1IY9HrfUZERLBu3Tpu3rxJrVq1shuexfEima2QkBAqVKiAlZVVhp8jZbb0MfW95YcffqBp06bUrl2btm3bUr9+fZPFkhJltrTkN7OVV6iUPQ9xdXXlpZdeYs6cOZQoUYJSpUoRGBiY6faZ6X3GxMTw5MkTSpQokd1w85Ts9B5T9j5DQ0PZtm0bDg4OfPjhh5QvXz6nQrQoXiSzZQzKbOljarPl7OzM/v37uXXrFrdv3zabInlltrRkZLZiY2Oz1Xk2V6SUWZ6eCJTZynMqVqxI+/btadCgAW5ubql+kf/55x8iIiKAhBF6fn5++Pr6MmXKFCZPnsy0adPS7X2uWbOG4cOH89NPP+X6+8hpSpUqRWhoaJbazp49m7t377J06VKOHj3KhAkT+PHHH1/YRAteLLNlDMps6WNqswVQpEgRunXrxogRI7h27RrmsJScMlta0jNbBoOB77//Po8jsgxUym5ikqZ+WL58OZMmTeLMmTMAHDhwgMePH+Pp6Ym1tTUGg4GDBw/i7Oyc5rns7OyoVq0a9+/fz5PYc5LatWtz+fLlLPVoa9WqxY0bN/jxxx8pUqRILkRneSizpY8yW/qY2mxBwj1g1apVPHjwgD59+tCjRw+TlwAos6UlPbNlZWVFYGAgV65cMcks7eaMMlsmRgjBsmXLOHbsGB06dEAIQXR0NKtWrcLT05NixYpRpEgRihYtihCCDz/8MM0REW+88QYXLlzIVrG5qahZs2aWk8SAgAA+/vhjlWilQJktfZTZ0scczFYSkydPxtvbOznRCgsL48SJEyaJRZktLRnVbL3++uscOXIkDyOyDFTKbmICAwM5c+YMrq6u1K1blx07drBo0SLatWtHsWLFUh3brVs3ChYsyMcff8y0adMoWbJkqv1CCAYPHszgwYOZOXNmXr6NbFO2bFnNwsGZ5ZVXXsmXNQLZQZktfZTZ0scczFYSBQoUwMHBIfnnZcuW5dks38+jzJaW583WuXPnOHDgAIULF6Z69erExcVhZ2eHwWAwatJuS0DVbBmHoRUAACAASURBVFkwb7/9Nu3bt+fkyZOsXr2agwcPMnjwYMqUKaM51t/fn6CgIPr27cvYsWMJDg5OtT82NpYTJ07w8OHDvAo/x7CyssryTW3//v0WafNyE2W29FFmSx9zMluTJk1i586dhIeHc/nyZWxtbTUdz7xCmS0tz5utU6dOMW/ePHx9falRowZly5alefPmXLp0yYRRmh8qZTcxSTVKFSpUoEmTJukuwNq5c2ecnZ1xcHCgT58+jBkzhunTp+Ps7Mzy5cvZunUrLVu2ZPDgwURGRqbqHVoCWbVTTZs2VWbrOZTZ0keZLX3MyWzZ2Njw008/MWbMGIQQbNy4kYULFwIwZcoUGjdunGexKLOlJaXZio+Pp0CBAtjZ2QHw2muv8dprrxEREcGnn35KjRo1TBlqjpLdARvKbJkJbdu2zXCle39/f27dugWAo6MjvXv35sMPP2T16tVcunSJQYMGUbFiRYoVK8a9e/fyIuwcxdHRMUuG6sCBA8psPYcyW/oos6WPOZktgOLFizN+/HhGjx7Nrl27GDp0KLt3787TRAuU2dIjpdkKCgoiMjJSs0qHo6OjKULLddRjxBeEzp07U7hwYaZNm8b8+fOJj4+nSJEiBAYG4unpmXxc0aJFuX37tgkjzRrVq1fn7t27mTo2JiaGBw8ecP78eapUqaLM1nMos6WPl5cXnTp1MnUYZsfQoUOpWbOmqcNIhYeHB/Xr18fT05O33nqLTZs25XkMymxpSWm2SpUqRcuWLXnvvffYtm1bquM8PT0t8nsot1DJloWwevVqBg0axNChQ3F2dubevXvs378fHx8fgoODcXJySj62RIkSRk2Wai7UqlUr3V/O27dv89tvv7F+/Xr27t3Lo0ePcHNzIyIiQpmt51BmSx9ltvQxN7P1PK+99hpVqlTJ81Fuymxpeb5mK2lN0jlz5qQ67q233uLYsWN5HV6uIYQgPDw8y+1Vym4htGnThpkzZ7JkyRKcnZ05deoUrq6uNG7cmBUrVqTqfRUvXpxTp06ZMNqsUb16debPn5/mfisrK1q1asWQIUNSbXd1dVVm6zmU2dJH1WzpY041W2kxZMgQJkyYwLVr16hUqVKeXFOZLS3Pj0aMiYlhz549/P3336mOc3V15dmzZ3kdXq7RpEkTpk+fzqxZsyhcuLDR7ZXZshCKFi3Kyy+/zKBBgwCoW7cumzdvJiYmRvPl4ejoyOPHj00RZrZwdHRMtxdZvHjx5Jq1lIwbN06ZredQZksfZbb0MXezlYSUEmtrawDWrVvHtWvXcvV6ymxped5s2dnZ4e7uztWrVzXHli1b1iK/i/Swt7enZ8+ejB49OkufCZVsWRDDhg2jcuXKyfandevW/Pzzz5pkSwhhsTeI9AyVvb09kZGRmu1TpkxRZus5lNnSR9Vs6WOONVvP88UXXzBp0iTKlCnD/PnzOXPmTK4bLmW2tOjNIN+hQwemTJmiObZHjx75aoLTokWL0rp1az799FOjRyeqZMuC8Pf3p02bNgQEBAAJX6iXL19OXk8xJZb6qMTGxibND3FUVJRm1Asos6WHMlv6KLOljyWYrZ07dzJ//nx8fX158OABH3zwQa5fU5ktLXozyNvb21O5cmU2bNiQanvVqlWzPFm1uVK2bFkqVarE1KlTjWqnki0LonPnzqxZs4b33nsveVvNmjWJjY3VTHBqqcmWu7u7rnaWUrJs2TLGjx+v2afMlhZltvRRZksfSzBbBw4coEePHvTv35+GDRtSvHjxXL+mMlta0lob8aWXXmLVqlU8fPiQcePGsX37dgCcnZ2JiorK6zBzFU9PT6Kioli2bFmm26hky4JYvnw5169fT7VMz5kzZ/jqq684dOhQqmPj4+OzPQmbKfD09NSty/L39+f999/Hzc1Ns0+ZLS3KbOmjzJY+lmC2AKytrdm4cSNeXl55cj1ltrSktTaiEAIvLy8+/fRTXnrpJXbt2gVAly5dOH78eF6Hmeu0bNmSw4cPs2fPnkwdr5ItC0BKyZ49e5IXanZ0dERKyfXr12nQoAGVK1fWDEl1cHAgNDTUFOFmCw8PD82C1IcOHaJOnTo0adJEt40yW1qU2dJHmS19LMFsAdy4cQM3Nzfi4+N58OBBrl9PmS0taZktSHjENnjwYFxcXHj69CkAjRs35saNG3kYYd7RuXNnFi5cqCsInkclW2ZOWFgY8+fPp3HjxlSvXp1OnToxceJEdu3axf79+3nnnXeAhA95yuTKUic2ff4x4s2bN3ny5AmDBw9Os40yW1qU2dJHmS19LMVs1ahRg3v37tG2bVvd0W85jTJbWtIyW88THR3Njh076NGjB4UKFcqX/45CCDw8PLh+/XqGx6qU3Yw5cuQIS5cupXHjxvj7++Pg4ICXlxdVqlThzp07FClSBHt7ewB69uzJn3/+yWuvvQYk1HJt2bKF2rVrm/ItGI0QItlShYeHExAQwNy5c9Nto8yWFmW29FHzbOljCfNsJfHw4UPeeOMNwsLCiI6OTr4H5gbKbGlJz2ylpFatWhw+fJiSJUvi5eXF2bNnqVu3bh5EmLeEh4fj6uqa4XHKbJkhkZGR+Pr64uvrS69evejRowfe3t6EhoYSHR0NQEBAQKrROFWrVk1lhFxcXLhy5Uqex54T2NjYEBcXx4oVK/jpp5+S59VJC2W2tCizpY8yW/pYitkCqFKlCjt27GDXrl0sW7aMkydP5tq1lNnSklmzVadOHZo1a4atrS0tWrSwmM+Xsbi7uzN9+nTGjRuX7nEq2TIzzp49S48ePbCzs2PevHk0adIk2dp07tyZcuXK8fDhQypUqICDg0Oqts//7OTklCm9aW4UKlQoeeRhsWLFMjxemS0tymzpo2q29LGUmi2AL7/8km7dunHlyhUePXpEZGQkZ8+eJT4+nuvXr3P69GkADAZDtq+lzJaWzJqtJGxtbTl48CDly5fPxahMR9WqVenTpw+vv/56usdlKtkSQtwQQpwRQpwUQhxN3OYshNgmhLic+FoscbuVEGKhEOKAEKJW4rZWQggphOiU4pwbhBCtsvoG8xsxMTH89NNPjBs3jmnTpjF8+HAKFiyY6hh/f39u3brFzp07GT58eIbnbNq0qVFDU82FFi1a0LFjRzw8PDJ1vDJbWpTZ0keZLX3M3WzFxsZy/Phx4uLiKFmyJEuWLOHGjRv4+vpy8eJFNm3ahJ+fHz4+Pskjlrt27Zrt6yqzpSWzZisJOzs7Vq5cSf369XMxKtOTUVJuTMreWkqZcjKn8cAOKeVkIcT4xJ8/BdoDh4BxwA/A24nH3wY+B9Sd7jmuXbvG0KFD6dq1K+vWrUvT0nTu3BlnZ2cuX75MkSJFNPufn+qhWLFi3Lx5M1dizk3at29v1PHKbGlRZksfVbOlj7nXbAkhGDlyJB4eHtjb2xMcHMy7775L06ZNmTVrFrVq1aJo0aJs2bKFcePGMXHiRN3VJoxFmS0tWTFbjx49okCBArkYlfmTnceIXYAFiX9fACR1I6wBQ+IfkeL4U0CYEKJdNq6Z77h27RqXL1/m3XffZfjw4ekmDUlmK61ffr15tUqUKMGFCxdyLF5zRJktLcps6aPMlj7mbrZsbGz4888/CQoKIigoiLCwME6dOsUnn3zCuHHjeO2113Bzc2PMmDHExsYyceJEevfuzYYNG9i6dWuWr6vMlhZjzZatrW2mn1LkZzKbsktgqxBCAnOklHMBVynlvcT994GkcvwtwGJgIPDec+f5DpgEbMtW1PmIM2fOUKpUKQ4fPkx8fHy6xeBJZiutOUv0ahSaNGnCsmXL+Prrr3MoYvNDmS0tymzpo8yWPuZutgCqVavGihUr+OijjxgwYACFCxemYcOGAMTFxXHp0iVOnjzJo0ePaNasGQ8fPiQoKIi2bdtm+ZrKbGkx1mzVrVv3hbdakHmz1VxK2QDoAAwXQrRMuVMmKBWZ+Pc4KWVvKWVTKeWZ547bAyCEaJ6ZixoMBs6cOZOvX48dO0ZgYCAeHh7ExcWle/y6devYuXMnQgj279+PwWDQvD7f7saNG9y9e5d9+/bpHp8fXseNG0dAQIDJ4zCn17Vr1zJt2jSTx2Fur7NmzcLf39/kcZjb64QJEzh79qzJ48jodfny5Zw/f5758+ezYMECtmzZgsFg4Pz589SoUYPatWvTp08f9u7dy+DBg5OH5a9atYpt27YZfX+eOnUqMTExJv+eMKfXe/fuMXfu3Ewf7+TkxJUrV0wed168pouU0qg/wNfAx8BFoHTittLAxXTatAI2JP69PbAZ2AC0SqeN3Lp1q6xfv36+fPXw8JCjRo2SQgi5bt26TLX7/fffpaenp/z0009l69at5bNnz1K9Dhs2TLfdF198IRs3bqw5Pr+87tmzR7Zq1crkcZjTa4sWLeSxY8dMHoe5vTZr1kxeunTJ5HGY2+tLL70k79+/b/I4jHmNiIiQpUuXlqtWrdLc96pWrSpHjhwp3dzc5Jo1a6Sbm5t87733ZL169Yy6T1erVk1u3rzZ5N8X5vS6du1aWaNGDZPHYY6vSe5J74/IaP08IYQjYCWlDE/8+zZgIuAFhKQokHeWUn6SxjlaAR9LKTsm/nwoMUEbKKXcnUYbuXu37q58wfz582nQoAF3797NcMhoEj///DPNmzcnLCyMCRMmaPaPGDGCnj17arZHRkayb98+vvvuu2zHbY60bNmSbdu25erkhpbGf//9h6+vL/PmzTN1KGbFwoULCQ8Pz9Ro3heJ0aNH4+3tTdOmTU0dilGEhoYycuRIBg8erHnct27dOpo1a8Zff/3FjBkzuHjxIgEBAVSsWJHKlStn6vyjRo3ip59+ws7OLjfCt0iuXbvGypUr+eQT3a/7F5pWrVohpRR6+zLzGNEV2CeEOAUcBv6RUm4GJgPthBCXgbaJP2eW74CyRhyf7+jUqROjRo0yqp6gc+fOuLm5pVmflFbi7ODgQFBQkEUuTJ0ZVM2WlvKqZksXNc+WPpY0z1ZKihYtyvjx41m9erVmX6dOnfD392f+/Pns2rWLdu3aERQUREBAQKbPr2q2tBhbs6VIIMNkS0p5TUpZN/FPLSnld4nbQ6SUXlLKqlLKtlLKR+mcY3eS1Ur82V9KKdKyWi8CkyZNYvz48Ub9Ivv7+3Pz5k2jRiMmUaFCBf7991+j47QE1GhELWo0oj5qNKI+5j4aMT1q1aqFl5cXe/fuTbXdysqK7t2788knn9ClSxcMBgOzZ8+maNGimT63Go2oxdjRiIoEVMpuAh4+fMjJkyeN/jLs3Lkz1tbWPHz4UHd/UrIVHx/Pjh07CAkJwc7Ojvj4eJycnFi9erXFPSbIDMpsaVFmSx81GlEfSxiNmB7dunXjwoULXL58GTc3N3bs2MGpU6f48ssvqVOnDosWLWLEiBFA+p3S51FmS4syW1lDLdeTx1y7do3Bgwfj6+uLELqPdtPE39+fq1evplk/kHQTuX79OhUrVmTWrFnMnDkTPz8/7O3tOXHiRL78olFmS4syW/oos6WPJZutJMaPH8+JEyfYtWsXQ4YMoVq1akCC+bpx4waHDx8GSF7WJzMos6XFEsyWlJL58+ezfft2Ll++TExMjKlDyrhA3lTktwL5sLAwNm7cSKlSpbhx4wY+Pj5Gn+PmzZsEBATw/fffa1YZj42NZfTo0fTq1YtNmzYxYsQIypQpk7w/Pj6eVq1a8c0339CmTZtsvx9z4tChQzRu3BgrK9V3SCI8PJybN2/i6elp6lDMijt37hAfH0+5cuVMHYpZcf78eUqXLm3UIzZzZsWKFQQHB1OrVi0gYRqh+fPnM3PmTBwdHZk5cyZXrlyhe/fu6Zqrc+fOUaNGDXVvSUFkZCT379+nUqVKpg5FFykl69atw8rKimHDhnHs2DFOnz7N06dPiY2NpWjRokavUpJZslsgr8gGcXFxbNy4EX9/f2xsbAgODs6ygv3rr7+wt7fXJFoAz549SzZeYWFhqRItAGtra1asWMH9+/ezdG1zRpktLcps6aPMlj75wWylpGbNmty4cYO///6b5cuXc/DgQdq0acOnn36KtbU1Y8eOpXnz5hkuZ6bMlhZzN1sXL17k6tWrREdHU7lyZXr16sW3337L9OnT8fPzo3Llyly+fDnP41IPo3OR48ePs3//foQQuLq60q5du2z1HB0cHPjqq6909z179iy5bimtKRBKly5N3759s3x9c0XVbGlRNVv6qJotfSy9Zut5ateuzRdffMHBgwfx8PDg3r177NmzhydPnnDs2DEaNWpEnTp1WL16dbrTQKiaLS3mWrO1bds27t69y71797Czs0tz1ZSBAwfy/vvvU7Vq1TyNT5mtXOD27dvMnDmTJUuW4OrqSr9+/ejZs2e2Eq3Tp08TGRmZZq1BVFQUtra2PH78GHd39yxfxxJRZkuLMlv6KLOlT34zWwDR0dH4+fkxc+ZM/vjjD27duoWbmxsnT57k+PHjVK5cmQcPHqR7DmW2tJib2YqPj2fZsmU0a9aMefPm8ccff7Bu3TpKlCihe7yNjQ3169fPdN1eTqFS9hwkMjKSpUuXsmnTJrp06UK/fv0oXrx4ts9rMBg4fvw406ZNw83NTfeYpMeI58+fz9ZaYJaIMltalNnSR5ktffKb2QIoWbIkZcqUoWPH5FmHMBgM3L17ly+//JJFixZleA5ltrSYk9mKiopi6dKlfPTRR9SpUwcgU5/jIUOG4OPjQ9OmTXFwcKBgwYK4uLjk6veIMls5gMFgYPv27QwbNozw8HDmzJnDgAEDciTRAti9ezdDhgxh7ty5afY+k8xWYGAg9erVy5HrWgrKbGlRZksfZbb0yY9mC8DR0TGVmbKysqJMmTIUKFAAHx+fDL9cldnSYi5mKyQkhCVLlvC///0vOdHKLPb29kyaNIlq1apRpEgRHj9+zJYtW3Ip0gRUym4k8fHx9OjRg8ePH6fa7uXlxfvvv0+TJk1y9HpRUVEEBwfzyiuv4OzsnGbWnlSz5ejoyPXr15OHPb8IKLOlRZktfZTZ0ic/mi2AV155hUuXLuHh4ZG87ebNmxw7doxWrVphY2ODwWBIc7ShMltazMFsXbt2jX///Zc5c+bg4OCQpXNUrlw5uV5PSsmwYcNyMkQNymwZibW1NWvWrGH9+vUMHTqUTz75hH79+uHj45PjiRbA5s2bk9egSq/3mWS22rVrx4wZM3I8DnNGmS0tymzpo8yWPvnVbLVu3ZoLFy6k2rZv3z5GjBiRPEH0s2fP0myvzJYWU5uto0ePcvXqVWbNmpXlROt5hBC53mFXyZaRxMfHs23bNjZv3kyNGjUIDg5mwIABmqkWcoLg4GAKFy5M+fLlgfTXL4uMjMTOzg57e3uKFi3KiRMncjwec0WZLS3KbOmj1kbUx1LXRswIJycnTbJUr149SpcuzdatW2nevHm6X9jKbGkxpdnaunUrjo6OfPvtt2nayIMHDxIVFWX0uXP7/1klW0Zw7tw5xo4dy8GDB4mNjaVAgQIMGDCAAgUK5Mr1Nm3alGpl9YzMVtI8W23atGH27Nm5EpM5osyWFmW29FFmS5/8arYAihQpQkxMDAaDAUiYUX716tVMnDiR6tWrp9tWmS0tpjBbSSMOX3nlFd5//33dJZfOnTvH+++/z6xZs9i1a5fR1yhSpEiWkrTMomaQzwQGg4E333yTwoULM2rUKF566aU0l8zJKS5fvkxsbCwffPBB8rb0ZnleuHAhQLIFCwgIoFWrVrRs2TJX4zQH1AzyWtQM8vqoGeT1yW8zyKfk6NGj+Pr6Urp06WSrKaXM1HJpagZ5LXk9g3xUVBQzZsygbNmyODg4cPfuXd566y369OkDwN27d5kyZQp2dna0a9cOg8HAxo0badq0Kbt372bcuHEZxhoZGcnGjRsJDQ3N1vxbagb5bCKEYPHixaxYsYLmzZvneqIlpWT//v289957qbZnZLbs7e3566+/+Pvvv3FxcWH+/PlGLbpqqSizpUWZLX2U2dInP5utRo0a4eXllWrEWmbXpVVmS0tem6179+7RtGlT2rdvz1tvvUW9evUoXrw4jx49YsKECUyePJl27drxxhtvYGtri729PcWKFSMsLIxu3boxYcIEnj59CkBoaCj79+/n119/Zdy4cXz44YcMGzaMdu3a4eHhwd27d3PtfSizZYb8+++/1K1bl9dffz3V9vR6n9OnTyciIoKKFSvSqVMnNm/ezJgxY5g7dy5vvPFGXoVuEpTZ0qLMlj7KbOmTn80WwOjRo+natWumk6wklNnSYuq1ER8+fMixY8eIjIykVatWGU6x9PjxY9asWYOjoyP29va4u7tTrlw5XFxckj8Pf/75J+7u7nz//ff83//9Hw0aNEhzUtT0UGbLgoiJieHKlSuaRAsyNlunT5/G29ubwoUL07NnTzp37syLkLAqs6VFmS19lNnSJz+bLYBly5bxySefsHz5cqPaKbOlxdSjEV1cXOjQoQPdu3fP1FyWxYoV4+2338bb25uuXbvSuHFjXF1dUyXePj4+tGjRgtmzZ7Nnzx4CAwNZt24dK1asYNWqVRw9epQnT55kK241zMLMWLlyJZ999pnuvvTmwjl58iRNmzbF2to6eVt0dHSuzx1iDqjRiFrUaER91Dxb+uTXebaSaNOmDQ0aNDB6omk1GlGLOcyzlRNER0ezadMmoqKiKFSoEM7OzpQoUYKXX36ZunXrMnr0aACePn3KsWPHOHDgAA8fPiQ6OpoCBQpQsWJFqlatSsGCBTN1PWW2zIiAgAA6dOiQpp5Nr/fp4+PD+++/n2rb999//0I8RlJmS4syW/oos6VPfjdbTk5ONG7c2OhHX8psaTG12coqa9euJSQkBICzZ8+yfPlyRo0axaxZs/jss8/o0KEDJUqUoHTp0jx69Ci5XaFChXj11Vf57LPP8PX15ddff+Xzzz+ncuXKHDx4kGXLlmXq+qpmy0y4efMmV69e5ZtvvknzmPxeV5FVVM2WFlWzpY+q2dInP99bnj59yueff063bt2MbqtqtrSYumYrK1y5coVHjx5x+vRpHB0dadiwIUOGDEl+lPjuu+9SsGBBBg0aRMOGDY0694YNGzh+/DivvvqqqtkydyIjI9m5cydfffVVusel1/v8/vvviY6Ozo3wzB5ltrQos6WPMlv65GezFRgYiLOzc5baKrOlxdLMlpSSffv2MXbsWCZPnsyECRN49913kxOtmzdvUqxYMTp37oy/vz/vvvsuK1euzHS5QceOHYmJieH27dvpHqfMlomRUrJgwQImT56Mq6trusem1fs0GAx0794dV1dX/Pz8UtVtvQgos6VFmS19lNnSJz+brc2bN3PlyhVq165tdFtltrRYmtnavXs3rVu31p1zMiQkBD8/P+rVq4eTkxOQ8J3833//ceLECfr27YuXl1eqNs+ePePy5cupPk+xsbEMGTKEhQsXKrNlrmzevJnBgwdnmGhB2r3PwMBAPDw8aNasGZ988skLMbdWSpTZ0qLMlj7KbOmTn83WtWvXcHFxyVJbZba0WJLZioyMJDg4WDfROnDgAK1bt+b27dvJiRYkzMFWu3Zt+vXrx44dOzTt7OzsGD16dKo5uWxtbfnhhx/SjUUlWybk/PnzuLi48Oqrr2bq+LTWL/vvv/9wc3OjfPnyVKlShW+//TanQzVr1GhELWo0oj5qbUR98uvaiAC3b982ehRiEmo0ohZLGo34zz//8Nlnn2lKbI4cOcKCBQto2LAhVapU0W1rbW2t+yjRysqKAQMG0K1bt1TndXNzSzcWlWyZiNDQUE6ePGnUF2Javc/Tp08nL9Pj4eGBg4MDP//8c47Fau4os6VFmS19lNnSJz+brWfPnmU5YVJmS4ulmK24uDgCAwOZNGlSqtVYoqKi8PX1pXfv3rzyyivpFsQnraf5PP3796dEiRJ8+umnma6VVim7CYiPj2fVqlX88ssvRs1onNZcOEFBQRQuXDj550aNGhEQEMDixYvp379/jsRsziizpUWZLX3UPFv65Od5trLz/63MlhZLMVs2NjaMHTuW27dvExERkbw9OjoaV1dXrKysMlwHMa1ky8bGhnLlylGvXj0++ugjgAwXsVZmywSsW7eOjz76iCJFihjVLq3ep17P69VXX+XMmTNs3LhRs+/EiRN8+umnRl3bnFFmS4syW/oos6VPfjZb2TFTymxpsRSzlcSpU6dSLVlnY2OTZhKVkpCQkDQnLA0JCcHe3p6KFSvSs2dPevbsycCBA9M9n0rZ85ijR49Sr1496tata3Rbvd6nlDLNROONN95gw4YNbNu2DQ8PD9q0aUNERAR+fn6pCgItHWW2tCizpY8yW/rkV7MlpVRmK4exFLOVRHh4OGXKlEn+2draOsNky2AwsGbNGmbPns2dO3cIDw+nRo0ayfu3bt1KrVq1jIpDma085P79+9y9excfH58stdfrfT58+JBChQql2aZjx4506dIFBwcH/vzzT/744w+8vb0zvcSAJaDMlhZltvRRZkuf/Gq2QkJC0r0/ZoQyW1osyWzFxcXh6OiYaltmzNaaNWsoWLAgH3/8Mb6+vkyfPj3V/sOHD1O5cmWjYlEpex4RExPDhg0b+O2337J8Dr3e59mzZ3F3d0+3nRCCMmXKJGf3ISEhWR4KbY4os6VFmS19lNnSJ7+arexMaArKbOlhSWbr4sWLmmkfMmO27O3tadKkSfJ36/Lly4mOjsbe3h4pJVFRUUbPvabMVhaRUmbquS/A1atXWbhwIRMnTqRAgQJZvuacOXM4d+4cERERREdHEx8fz+nTp42eoDE0NBQpJYsWLUpeK8qSUWZLizJb+iizpU9+NVvXr1+nRIkSWW6vzJYWSzJb58+fp02bNqm2WVlZZfjd3bFjx1QSo1q1auzduxdIWPqnZMmSRseiUvYssm3bNq5du0bHjh1TPQ9O6I6bEAAAIABJREFUyf3799m2bRsNGzbkzz//zHYPaejQoTx58oQ333wTT09PDAYDVlZW1KtXL9PnCAkJ4eDBg6xevZoWLVrQuXPnbMVkDiizpUWZLX2U2dInv5qt69evG/24JyXKbGmxJLMlpdR9jGzMLACQMKXSjh07cHJy4ueff6ZXr15Gx6LMVhaIjY3l8ePHLFiwgIsXL7J9+/ZUs7aHhYWxfPlyrl69yowZM/jggw9y5Bd2zpw5xMfHU6FCBXr27Im3tzc9e/bM8IMTHh7Oli1bWLp0KZcvX+bLL78kJCSEtWvX5otCeWW2tCizpY8yW/rkV7P14MGDbC1BpMyWFksyW/Hx8ZoVVaSU2NvbG3Uee3t7rl+/zqpVq/Dx8Uk11VJmUSl7FtixYwcffPABdnZ2TJw4kZ07dzJ//nw6duzIgQMHsLGx4euvv87xuqi+ffsye/ZsGjVqlOk2BoOB5cuX8+WXX1KtWrUcjcdcUGZLizJb+iizpU9+NVtxcXFGW4yUKLOlxZLMVvHixbl+/XqqdRyvXr2apWTpgw8+SHd/aGhouvuV2TKCb7/9luXLl/PkyZNUi1C2adMGX19fTp8+zfDhw/nxxx9zPNF69uwZ/fv3p3r16jx58iS5t3X//n3++OMP9u3bp7sm4q5du/jggw/ybaIFymzpocyWPsps6ZNfzVa5cuV4+PBhltsrs6XFksxWjRo12LVrV6ptGzdupE6dOpk+R2xsLBMmTODBgwfJ25KK5IODg7l58yYTJkxgyJAh6Z5HmOuixUIIuXv3blOHkYp//vmHb775hvj4+GwNJ84K8fHx9O/fnzZt2nD27FlCQkJwc3PDysqKCRMmsG/fPpYsWUKTJk2S5wOJiIhg8+bNzJw5M09jzWsOHTpE48aNjR4dkp8JDw/n5s2beHp6mjoUs+LOnTvEx8cbPagkv3P+/HlKly6drUdu5sjFixdZvHgxbdu2zVL7c+fOUaNGDXVvSUFkZCT3799PZYvMFSkl69atS9XxHDVqFN27d8/0OdasWUPFihVZv349Q4cOZdu2bZQqVQobGxuKFStG0aJFcXJyomDBgrRv3x4ppa5KVZ+gTBAXF8eKFSto2bIlBQsWzPNECxKGq547d46zZ8/i5OTErFmz6N69Oz/88AMFChSgbdu2zJs3Dzs7OxYuXMi9e/fYsGED48ePz/NY8xpltrQos6WPMlv65FezVa1aNYKCgrLcXpktLZZktoQQqdYujI+PJyYmJtPtQ0NDmTFjBsHBwZQvX56LFy8ydepUJk6cyFdffcXIkSMZMGAA7du35+TJk+meSz2MzoCgoCD8/f35/PPPqV69ukljqVOnDm3btiUkJIRhw4ZpEikrKyt8fHzo3bs306dPp0qVKhmuRJ4fUDVbWlTNlj6qZkuf/FqzJYSgYMGCSCmzVLulara0WFLNFkCBAgUICgqiZMmS7NmzxyirvWnTJqZOnUpQUBA//vgjxYsXT7XfYDCwePFiAgICMrSnymylQ3BwMFu3bmX27NkmT7QgYQLTW7duYWtri4uLC++88w5Pnz7VHFewYEE+++wzxo4da4Io8x5ltrQos6WPMlv65FezBdCwYUOuXbuWpbbKbGmxJLMFUL16dWbPns2IESPYtGkTL730Uqba3blzh5CQEN544w0mTZqkSbQCAgKSv4MHDhyYodhQNVvpcPbsWcqVK0fHjh1NGkcS48aNIzw8nPr169OhQwdcXV2JjY01yWNNc0LVbGlRNVv6qJotffJrzRYkPJ344Ycf6NKli9FtVc2WFkuq2YKEMqC9e/fSvHlzo56AxMXFcePGDQIDA3nw4AHW1tbY2dlRrFgx7t69i7u7O82bN09lTFu1aqVqtrJCTEyMZl0lUxIbG0uPHj3w9vZm+fLlvP322/j4+HD48GFTh2ZSlNnSosyWPsps6ZOfzVahQoVS1e0YgzJbWizNbNnY2NC6dWujS01sbGyoUqUKrVu3pnfv3vTs2ZMuXbrQsGFDunXrRosWLYx6NG1RyVZcXFyWdXBWiImJwcHBIc+ulxFDhw5l1apVfPTRRxQsWBAXFxcGDRrEtm3bTB2aSVE1W1pUzZY+Xl5edOrUydRhmB1Dhw6lZs2apg4jV9izZw9Vq1bNUltVs6XF0mq2chpHR0fs7OyMbpfpZEsIYS2EOCGE2JD4c0UhxCEhxBUhxHIhhF3i9kJCCH8hxE4hhFvitsFCCIMQok6K8/0nhKiQmWsfOXKEv/76i02bNnHy5Elu3LhhxFvMOpGRkRQsWDBPrpUZJk2aREREBAMHDiQ6OpojR45QuHDhVPN/vIgos6VFmS19lNnSJz+brd27d+Ph4ZGltspsabE0s2UuGJOyfwicB4ok/vw/wFdKuUwIMRt4B5gF9AfmALeAUUDSkLnbwOeAtzEBbtu2jVq1ajF27FisrKyIiYlh2LBhVKhQwZjTGIXBYGDz5s3Y29ubVW8vJiYGb++Ef75GjRolzyQfGxtLfHw81tbWpgzPZCizpUWZLX3UaER98utoREhYPs3Y5VmSUGZLy4tutrJKpsyWEKIM8Cbwe+LPAmgDJKW3C4CuiX+3BgyJf1I+0NwA1BJCZHpY3/79+6lcuTJ9+vRJLlC0s7OjbNmy2ZoVOD1u3brFn3/+Sc+ePZk4caLZfImHhoZy9epV7ty5o9lXtmzZDOf4yM8os6VFmS19LNls/XvvXyYenIj3Bm86rulI27/b8uryV5lyZApRcVHZOnd+NVtPnz7NVnG7MltalNnKGplN2acDnwBJCwoVB0KllEmfwtuAe+LflwBL+X/snXd4FNX6xz9nZnfTA+k9hBZCQuhdBKQ3FUHEawNUEAXEn1fEawX1ylVRUVARywVEkKogCIJAKNIhAgk1CQkJpFCSkLpl5vz+2BDlJiAgklX38zx5TrLT3pndzJ75nu95X3AFHvzVPnTgLeB5YNjVHNTV1ZURI0Zc8trZs2c5ceIEXbp0ucrQrw6bzcaqVasICAjg888/d5hO1kWSkpLo0aNHtbOo4uLi+PHHH2nVqlUNRFbzOJWtqjiVrepxWGVr+3ZK167HvVd36NABNv4A2/fCbbcBcHDJdGa47yS1cSDPrIH4xFzO1AkkL8yd/Lkf8EKnZdz/8AxapJShbN4CXbsCoK9dj3Jxn9u3o2/ciFKxT7lxI3TtiujYkVFt2lJ39Wp7LB061MQV+EPw8PC4piSW/4tT2aqKU9m6TqSUV/wBBgAfVfzeFbtC5Q+k/GqdCCDpCvsYDszA3rk7CtQFkoCoK2wjzWaz7N+/f2VbXl4uIyIi5LJly2T79u3lunXrKttWrVrJFStWVHn9atqPP/5YhoWFycTExEuO50htbGysbNq0qYyOjq72PKKiohwizppoO3XqJPv06VPjcThS27lzZzl8+PAaj8PR2mbNmslp06bVeBy97rxbvvvDYdlz4D3SnJAgy1SDtArF3n4wTeoGpBRIq4K0mYzSKpBmF1VOi4+TmqJIHaQmhLQJRWpCkVZVka89ECItqiI1gbQqQmoGo9SFInU3N/lBfLy0GI3SJoTUTSZpVYXUFGS5osjSGTNkqclF2hRFlqmqLE9IkD2H3iu3ncmXfe+40yHet9/TTpw4UbZp0+aavxfWrVsnvby85KpVq67re+Wv2s6aNUsGBATUeByO2Nq7VNX3adRJkyZdsTM2efLkEcA9kydPfhK4F4gBooAmkydPfnfSpEn65MmT44BWkyZNmneZfTQHQqWUKydPnmwBhgCRwPxJkyZVWyp78uTJkyZNmkRQUBD169cnKCiI1atXExQURPPmzfH19SU0NBQvLy92796Np6cn5eXllJaWkpmZia7rAISFhVGvXr3K9X/dBgcHk5ycTGhoKI8++igdOnS45HiO1Lq5uRETE4OUktatW1c5n9TUVO69916Cg4MdIt6b2fbr14+WLVvSoEEDh4jHEVo/Pz/69OlDfHy8Q8TjKK2Pjw/du3cnOjq6RuP48aw3avJiugYVk/v1RqJTj2OQ9nvW0SNHCczPBwmKBKHrqBIUFBqEheOenW33ZwiBkBIFCQjiM9xxLyu2e0MqthNI0HTqBgfhmZ2DKiXoOoqUKBIUIVDMZgypqai6jhQK76meLB8xjoV5BXi0u4UH64QSEuwY79/1tG3atCElJYUWLVpU+z1wpdbV1ZXY2Fj8/Pyuabu/chsVFUWtWrWIjo52iHgcqZ07dy6TJk2aXF2f5pqSmgohugLPSCkHCCEWA0t/ZZA/IKX86DLbDQdaSynHVsxaPIR9SLKdlDL9MtvIX8eWmJjIl19+eUliuqSkJPbt28dzzz1XZWpvWVkZycnJ/PTTT+zbt4+HH364yjGysrIoKSlh1KhRV30NapJ77rkHf3//SpP8r0lISODuu+8mLi6uBiKrWTp37sy6deuu2wT7VyQpKYn33nuPzz//vKZDcSjmzp1LUVERY8aMqdE4Hpm9m9dO3EOoOA+ZNuTcUtBAqEAfF+QPFoSG3fVqMIJNA5MJpk2Dp54CiwVUFYRA2mxowsDmDk/QedsMFF0DgwGQCJuGVI1k3Hof4VvmoeoaomKZtFkRLq6Iae9jGTMGg66jmUzMmbsQQ8eOFNk0/p2WzbjIQF6o/+ct+yWl5PHHH+cf//jHNW/75JNPMnXq1Oua6v9XJS0tjSVLlvDss8/WdCgOx5WSmv6eweiJwNdCiNeBROCq7upSSosQ4gPg/Ws52Mcff4yrqytvvPEGY8eOZcWKFXTs2JHPPvus2sRibm5utG7dmp9//vmy/q7U1FQGDhxY7TJHxNPT87L1l86dO/e3qINYHU7PVlWcnq3qcQTP1uHsC3geW0awqQA6jIX7uyI83kFuOYytZ3sMjwxHdNoKSxZDl3io3w2+XYE+eBDKxQfDpUtJdcumPD2bgMJiEnv9kyNB7UC/QP2M/aTWvRWAtmV70eO7oEXeSpKLhtexraTHNSeulcS0cg3eI17HOGoU57KzMf74HSd6tmT44H4YFPtXQ0aZhekn87jVx4vOvl6XOSPHRghBSEgIRUVFeHld2zk4PVtVcXq2rg+HLtfzv7FNmTKFvXv3EhYWxosvvkhAQMAV9yGlZOTIkTz44IPVLt+9ezdWq5Vnn33W4dMmaJpG8+bNefzxx6tVrxYvXsyMGTNqILKax6lsVcWpbFWPIyhbixd9yZBDY8mp1ZzgRxdC0glk9+5IsxmpKuy7oxNtlieABiigKwbQNCwGwfzRbRk+cx/CZkXqEoQACZrBwKbWfeiyexWqpiMRoKgIqaMZVRaMbsu9M3disGroioJQdBRNgosbYto0yseOxahZsRgEX894lBEjZwFQqum023GIW328+Ci2To1ds9/L0aNHmTt3Lr169bqm7ZzKVlWcytbl+UuU65kxYwZz587lscce4/333//NjhbA3r178ff3v+zyNm3aEBwczCOPPEJKSsqNDPeGs3HjRvr27VvtbMSMjAyaNWtWA1E5Bk5lqypOZat6HCGD/KC0lwB4PO8uXlp/Fuv6DWA2o+g6qs1Gm8Qt9o6WxD60aLOh6BKDTSdu7UGExYLQJQIqPVuqptEw8wSqpqMAAomi21CkjmrV6LjyEIrVZvdqaRqKTSJ0wGImf8HXGDQNVQejTVL4wwqsmj2ViruqEO3uysmy6yt34yg0atSI3NxcrlVccCpbVXEqW9fHn6KzdebMGXJzc0lMTKRnz55XvZ27uzt79uy5Yl6devXq8cADD/D++++zfPnyGxHuH8KqVatITU3l5MmTVZbt2bPnb/3hd+bZqoozz1b1OEKeLfX2aQAsdn0d793v89wZb6RBBQHlqol3ou+lTDWhKwqYTAiTC1JVUVxcCX/4/9BUgcTeF5NCQRcKmmIktV43rAYTulDQFQVNUdCEQFOMqLYgFLtdHgDdINAV0IVKpuYNioouFFCM+Lt24Z/r/w+rbmXvmnR8SnUyy68/fYKjMHjwYLZv335N2zjzbFXFmWfr+vhTDSNeC0lJSbz99tu0bt2asLAw/Pz8fnObuXPn8vHHHzucZGy1Whk3bhwdO3bE19e3iu9g0aJFfPjhhzUUXc2zc+dO2rRp87uSF/7VKCoqIiMjgyZNmtR0KA7FqVOn0DStWoX4plKUC6ufhUPfcoxI1mbfymOpXzE05CXenPoYoYd/xmP71sp8WSQk2H/v0IHU+zvhtfYYOX7unPOCLH8rgVonjBZ3Stw30cCjHwdLFazGI1wQKQTI/yM9wkLUkTn4F2scjlUILCoi3r2U883cSXTtRP28o9Q5kcn5pu6cj/GkILU56cWtqJ0Rz+ZYVzbFu5PeuSmu6p/7f+zxxx9n4MCBuLq6XtX6hw4dIiYmxnlv+RWlpaXk5ORQr169mg7F4fhLDCNeKzNmzODBBx+kadOmV9XRArj11luZNWvWHxzZtbNmzRri4uJYsWJFpbJVUFDAkiVLWLBgAV0v3pD/pjiVrao4la3qcQRlCwCvILhnDtw7nwhTMf2sP5Ch+2PARpS/Bx4G4FwKJH0D2z+i5OevYMmjMCUCW1I6ACktrOS3dkdTbuGTwNaUh6QTn6xzJn01Ui9FlJajCp1TtY/gam2O9HyEEtcobLINOU2CSTwnwKwT3Wg7Pr6nuYCNfF1HL/ckOO4nGob9iIUD/NzQirdiQ61mItJvIXWJteT3Zbe/kTz33HPX9P47la2qOJWt6+Mvq2y9/fbbLFy4kPbt2+Pu7o6LiwudOnX6zSeaefPm8cEHHzhMAert27cza9YsHnroITIzMyuVrf3799OwYUP69u1b0yHWOE5lqypOZat6HEbZ+hW2TydieOItNE1gUw2UvjUVn4lPgU2CApoAqYFiVOHhrjBrA+gSq6oCCgZdQxqMbO3+MJ3WzkTVJZqiAhJFSiwG+HhwLGOWHMOgaejCPmQpdBs2gyA1LobopOMomoamGFl++5tot2/FP2ofX5ifYL1rd0Iz3mPt0Hfwd7u8B7Y6Nnz1AjLka/SizvS8879/wNW7dv7zn/+QlpaG0WikoKCAhx9++LK+LKeyVRWnsnV5/rTKVmJiImPGjOH5559n/fr1mM1Xb9KcMGECU6ZMoVGjRvTp04ehQ4eyYcMGFixYcMUnla5duzrMkNy8efNYvHgxw4YNQ1GUS5QtV1fXa7oef2WcylZVnMpW9TiMsvUrDDmuSA1UKTFoNnL+uwDdBkjQdYGiCQwSdJvO7vXnkPY+GCZNw6RZUaUONiveR7YipD01l6JrKBXJS0026HwgG0WzokgdRbeiaFZUKTHawOOsu900L3VU3UrI6Z/5tLSchJRnWO/anR76Gp61mln47iTeWDYJKSUWi4XNmzezbds2fv75Z1JTUzl//nyVtBpFZfkAKF6bWbfo6UuWlZXlUlycc3Mu8q+YOHEin3zyCQ899BD16tUjKSnpsus6la2qOJWt68Ohp1ksXryYQYMGoWkaBw4cYOnSpRgMBiIjI+nfvz8xMTHV5ti6SM+ePS8x1L/55pssX76co0ePXjb5Z3h4OJs2bbqunCzXS35+PsePH6dt27YA6LrO66+/jtFovCSJ6x133IGvry9g72wVFxfflPgcHedsxKo4ZyNWjyPk2apCjz7w+htIqw2rauB0ODQ+LEEDqUhsGFB0HV1RiIzOQT0hkRcTngpAB1XVMcRa0bMU0HRsigoIVF1DU1WORLYn/thaZEXHQQqBhsRmUEiteyuh2QeRugXNIMgKj6PZmbZ8casXtcpO0bd0HQH10wmoD2azG4vnb8NaEEHGYQXN3RvNzdOeggJ7TqtaXt74+Pri6+dLZLOHOZUscItcg+K/nITvzZRbTyINpzC5FWKzGQkwjKdlt8dv2uW++J1x5swZtm3bhp+fH82bN692XedsxKo4ZyNeHw6tbPXu3RuDwYCLiwstW7Zk6NChDB48mHr16rFgwQKefvrp397J/9C9e/ffrG7fs2dP3n//mnKuXjf//Oc/6dKlCydOnADsEu2YMWMIDQ2lY8eOl6z7v8qWs7Nlx6lsVcWpbFWPIypbdOiAePVZRF2VhJ7taNLmNEd71eNEg3pofT0w9VVR6iqY+qh4tpL8NLQJyU0C2Ny9Makdm2NuGElu1zgiWtow9jGg1lVw6avi0lfBUFfg0kfh/nZbMfVRMQe4kBsUTF5gEHm+Hnx7Sxjb/K1s7TCa7NAYjL0MRLdewZJOboQVF7B83/PUKrcibUZKcwOhIARv/yyCGycQ0/4QXmeS8Es5TONsC62LAmlqicQv35XitHMk7T3A8lUr2ZMegLnE/pCoua4Bl1TMJbUpyeqAZnMhzzqdpG1zb/plDwsLo3v37vTv3/+y6ziVrao4la3r40/XZc/Pz8dkMtGoUSNKSkquefuUlJTf/GIODAxkw4YNFBQUULt27esN9TfJysri2LFj3HPPPQwdOpSsrCyef/55Bg4cWK2p36lsVY9T2aqKU9mqHodUtrZvR5/0HsKs0efkVlBd8F2rITQdTqhIIRE2G3oWeFFAhzV5CF3A0UJ73iibjYAUhU3mJ+m87UMUXYcMewdB6DrypE6xbRBeP36Fi2YlSOZU1lUcsqmU3MAlBJw5iyJ19FyFXY3qYDYY+TTlWWJs6SxPHkVZihseHvnU8yrnVEEYwUYNn6gdxN57ntKD7Sk5VEpO3h5KbYU0aduDvo36QaEk62g6J615nE1vSWjseoSQ5J+tS3j5YDoNG87J1M0cOjGOrAtvYtjnQUzLwTftstevX58TJ05wyy23XHYdp7JVFaeydX04tLJ1ESklBw8eZN68eSQlJZGYmMj3339/zckJZ8+ezaxZs1i3bh1FRUVXXDc6Opr9+/f/nrB/k0WLFjF06FCOHz/O7NmzefHFF3nggQcuO3vy18pWZmYmUVFRf2h8fxacylZVnMpW9TikspWQABYLoiKJKYdsGDQbBqkjbDak1YaQ2IcAD9lQNInQdfs2ViuKlAgpiTm6qiKRqQRNQ2gaQkqEpqEf3Wl/rcLTJaQ9KSpS4ivzUaU9AarQdCL32ScHvRFi76zfpXxPf7GGZiUnWF5s5VP3vXxICvszHmWM4RMutEzF/Y5teA1II3hgAeeCl5FQMJ5N+v+R3uhjtLhluNfOofBMFACBoUew1Ps3SxbdCcIdd1MnjKZyTuRO4lTKtpt22b29vbntttvYs2fPZddxKltVcSpb14fDd9kPHDhAcnIyvXr14vHHH/9dZXW2bt1KWVkZb7311m8+rQQGBpKWlnbZuoo3gpSUFFq1asXw4cPZt28fw4cPv6IH7dfK1r59+/j000//sNj+TDiVrao4la3qcUhlq2tXhIsJW7kZjCqGofci3v0aLFYURWLFgKrrWA2gxhrgpIauYU96SoV6ZVD4uWEdep1NBeyJTkFHkTqoktpxaZCp2zPSK6Ch2LdDcF74EEQuIJFCUOgXj0+xxqaQxvS3LWZ04QwKst3INnshZQbe7m6c8XbhK+9GFCnePMsHPCf+TUOXVBTDxRnkEoRECB2hVPyORLcZUQz2ByPfkEOknPqHvfcHmFxK2bf/XwSEr8Lk6nlTLv19993H2LFjadiwIbVq1aqy3KlsVcWpbF0fDv8pOnz4MDNnzrwh+3rppZeIiopi1MVirlcgICCAn3766YYctzry8vIqa/kpikLr1q1/c5sVK1bQrVs3XF1dad++vXM6cgUTJkxw1kb8Hy4qW87aiJeyfv36Gq+NWIUOHRDTpnFw2ucsjGzDmehhvPBSC6LWLUHxTwNzMYnJERhjNbI72OjDKZRDNpTYitv3IRuWxi60aXUEg4+BsmR3doU/iKtaQHz6SvIbu3KwyVD8XL5DOXwac6wRq1sdYrZc4PP2sLmlF098E8wtRzOwGF2wmDyJOGMj39OFvZGBjORVXKLL6LV5OXHHz2Jx7U9q466V4Qup86bXv4g/uo8OezbibStG8Tag59lnSxsDVWo39cBFqASFBRHsXRf9pOCk5yyMrvYRhmZN/8v+rctw9/mODd+Oo8+9Ny9NxGuvvcbTTz/NsGHDqjzszpw501kb8X+4qGw5ayNeGw7d2Vq0aBHnzp27YfurU6cOTz31FBs2bKBZs2ZXTHbq4uJyXZ6wq2XRokWVsw+vlovK1qZNmxg9evQfFNmfD6eyVRWnslU9Dqlsbd+OfOopmpktNEk7wBSTgZDVH6NrNjRFwYiNVvphzJlGdjEI7Yd0VE3HllExFKhLlJNlnBbgu9qCq2amU8YMQEVIjcBMI+FmaLQ5HaHZ0I6BKg6iapKJGQKLKZSee09jtEoQCq0SP2FT83c4UNf+8HLHD/PZ2aoX3/W899K4pQ2EASkUPGURyY2akxzdnBZpicQn/Yw/GQBYTaG4izKsFhd+PqBTUJpL+w4d6d07kYKC3Xh7N0VVXalV+yCFOliNaTf18vv4+DB06FA2bdpUJUG0U9mqilPZuj4cWhoZNmwYH3300Q3d58svv8ypU6d+czhS1/U/1Ad0+PBhQkNDr2mbi56tzp07M3v27D8msD8hTs9WVZyerepxVM+WNJtRdA1htfDEuZ9x0SwYpI6iaUgNhAQXzUrPQ1tQNB0kKJrdY6VKe5oIl4NuCE1W5NmqyKdVkTvL69hmVE3HIMGggWK1F6KWQtJzXyFGq8Q+AqizL/IEjyz9D/VO2x90V/S+j1z/apKZCnsnRNVtlEk3vM1FNM0/wp4GbfjvwJFkt2pEVlhDyrpKCuplERjzI/FtlnPLrfNBPM3q1c35aesTfP9dN75b3oFC/V3MZV74u3a/aZf+In369KGkpIS8vLxLXnd6tqri9GxdHw7dZW/Tps0VPUzXQ3l5OW5ubr85y/Dnn3+mT58+17TvsrIyCgsLsVgslT9Wq7Xyd7PZjNlsprS09Jq9Z0VFRTRt2pTIyEi8vLzIyMhASnnDr8+fEaeyVRVGOx/yAAAgAElEQVSnslU9Dqlsde2KcHFBt1iwKQamejXhVcN2jJoVRdg9UFIHqUJU3DlkJuga9krUQqBLsKoGZjfqySun0jFpViQVT9LSXrJauljQFRV00BWBLnRUXWI1CNa19qbtsRJ0m/3vH9p64llYzONL5jHhyfEEnDvPGT/fynBjzqXRfOtGXMyl2FwC2dyhC5mBweS7+pDv6lO53rw2D15ymtEFu3kofxtBhgvoBjNCtSEUK6gWhGLDktuBTr1ex9sv6g++4NUzadIkRo8ezYgRIyotGk5lqypOZev6cOhP0R/RkXjvvffo3bv3b6536NAhxo8ff9X7XblyJUuWLCE4OBhVVVFVFYPBUPlz8TWj0YjBYKBfv37XFHdiYiLffPMN6enpxMXFER8f7+xoVeD0bFXF6dmqHof1bK1fj0hIQO3UmThjOE9ENSL68B4ONmhOn6ALdD2xHf/aB3EJKiY/2MjuXBd213OnsdaYu0K6cLRBc4w2d9529SXyXCYr6sYRnZ9Ov8wlRBa4khxmJTEilNBCf3Ii76dcQETWHvaFHWN/oyB+ngo++8s40NiddJMvnyzP5vt/xgDQe1cKjXwKeam9PUH0Eb96HLmzHg1zz6KbLyBMlw6QeNpKGXBuNxEBW3ARZmS5gQVlLTjmM4Ap7k0ZdkHBfHgNQQ3rkuPtzqReAzCpNe+Jcnd354knnuDbb7+tLIPm9GxVxenZuj4curN1o8nIyMBisVw2M3xhYSF5eXl4e3sTHR19VQZ0s9nMK6+8goeHByNGjLih8UopOX78OCkpKWRkZNC3b1/y8/ORUjJhwoQbeqw/M05lqypOZat6HFLZAqw2HaFLjAaFB9vX4Z7iaC5sOIdPv7YgdJSPvkOxFcPDP1Br5gv0TttGvbII6th2cDa/Gc3Ca9O8Y0f2f/0lpVaVgel7WR7ZkqyctoSJA4ScF3hZ3fHKP0Xt4tkQ3AKPU5tJ946k+/EoShQ/fPiOCJNkfIgHOdPOUwu7d8q1+17qGX9koF7It8ovisbxIH98dRV/kcNAuYkmHKAhxzCoNgj85dykC0w8bGGiJZEo0wQ+8nfBtV0vRHoJelI5LoZVvNLzrpt8xaunffv2rFu3joyMDOrUqeNUtqrBqWxdH3/ZQtTVsWTJEjIzM2nZsmXla3l5eWRlZdG8eXOWL19ObGwsmzdvZvr06fj4+Fxhb5CcnMybb75Jv379CAkJueZ4zGYzP/74I40bN65S1PP48eMsX76c4OBgXn75ZXbt2sXEiRPp06cPM2bMuOZj/ZXp3LmzU9n6H5KSkpzKVjXMnTvX8ZSt7dvRunVDmi1oBiM7xr1E5+kvI2w6GFQ0qSM1iWZU+G5sdwZNW29P26BWZE3QQVONbGk3ik47ZqJq9hQOAomQOlIIcgKDCM7LRUiJxJ42QpESm6Jy/IEYYr5ORlgl0ihY93YTsuICieIE/2QG3hQyig8xYOP/xMcA9JErWU9P3ChjNNOpt8oDt5JMDPecvOTUTm5pSUj7ZLRSX37MCOaxs4JpvreSUieEJN+GKFYbplN57PpHNwJdXW/6pa8Oq9XK2LFjue+++3jyySedytb/kJaW5lS2LsOVClH/rbrsgwcP5tFHHyUuLo7S0lJ2797N119/Te/evTl+/DjFxcWMGjXqN1NDSCn55JNPOHz4MMOGDav2ycdisbBz506sViuKoqCqKr6+vpfUZNyyZQvJycmcOnUKIQQjR44E7Arbjh078PX1pUmTJkRERBAREYHZbL7m4ce/A05lqypOZat6HFLZSkhAsVoRFQWlxbIlSJu0Jzm1aihUJCK1Seqs3YWUut2PVXEaAlA1K42ylqFq9vUvPqcKQEOQ4hZMILmVN3xV1xGAUdfw3XUGYZUoOuhWyYn98UxuMoHOcgPhZLJXtGUC04mSabjKMsqFG2vEAO6RX7GdTrwlXuLFbtNorGZjxYCCjooOQOSt+wAw1sqmh6UnfvIbninI4DNrHYLcG3EkvAWnourTfNthOnt68EhEAD2Ca9eoRcJoNDo9W1fAqWxdH3+rT9GRI0c4efIkmzZt4syZM+zevZtPP/2Uc+fOMWTIEK5GSTt37hzPP/88TZs2ZdCgQVWWWywW1q9fT3FxMcOGDcPX1xebzYbFYmHq1KmXdLZKSkqYOHEiixcv5vTp0yxbtgxFUSgrK2Pq1KmsWbOm0jsAsGvXLurWrUuHDh1uzAX5i+D0bFXF6dmqHof0bHXtijCZwGLBYDLReOzDWP41AYPNimJQURQFabOhqQbSGvWm6dHFSE1WmOCFvWelCkLjz6NnCqQu0RSBQNgVMCEw+dtQMiVSBx1RmT1eAvkugQQa89CtEs2gUtLM3tHYLLpVhqhIjXRxqfq+SNxPD7mGTOrwuvtTNJP7OEQcHpRynz6HqOw8wsKOVq6f53OYz86OQlNycC1TqFtupf65/Zz1SOVghB8/aU1JKCmlSWYuD4cHMDDIB3e1ZibMZ2VlkZiYyBdffME777zjVLZ+hdOzdX049DDi0aNHiY6OvmH7lFLy2GOPER4ezosvvsixY8eIiYm56u2PHj3KG2+8wZAhQyp9X0VFRaxatYpBgwZVdrI6duzIjz/+iJeXF02aNOGhhx7i2LFjLFu27JI8LrquM3PmTD755JOrqsF4+PBhQkJC/tB6jX9Gdu7cSZs2bZxJXn9FUVERGRkZNGnSpKZDcShOnTqFpmlERkbWdCiXMmsWculSxODBMGoUZ9+bTupHc1jeoD1djIdpkHSMT6PvZGnzngzev4ZBGRtY09aMLqz03HOBda28iIhoz73rfkRP0in29EIToJh1TvoG4WWw4ltQgFdREUVeXlgUBb/CQgq8a7EsvjsdIn8mMjGN1bV7srFdKw52bccZgze38w1pNOAQ8VjELw8zt5Qn85Nr3CWnECxP05REjhFDuqhPU5nIcD6tyE4PmqaSltaanOyGRNvCaGsN4d9Ry1GLILw0nGKjjcxAP4407EaulLT2dsdNVeju683oyEBuJhaLhY0bN/Lll1/i5uZGWFgYrVu3vqzf9+9EaWkpOTk5VawvTq48jOjQna2XXnqJV199taZDqWTt2rUcPXqUpk2bVr62cuVKCgsLyczMJDo6mqioKHJycioVqX379hEUFERWVhYxMTGXdJQ2btyIoii88MILV5UK4qmnnmLo0KFOZet/cHq2quL0bFWPo3q29G7dkWYzwtUFfeo7KE+NQ9g0bIqKjoKqa1gNCivuCOWub0+j6DqaogDSvswoWHZ7DPcvP4ymCTShoCBRK5St/Q3r0izlBEKX2L8K7IWopRDkBgYRcO4Miq6hGVQW3v0UmWGhfNKzKxElGiOL3iQgLJm5PMxG0fOKp9JJJtCcfWQQxTr6oKEykKUMYDm6xYDJVM7586Hk7R9Il7I2HCw/wIf1a1PLtJI2Z2JwcdtJl5ILrBr0OZ/nlwPgpSrs6xiHl+H6S7VdL507d2bt2rUcP36cb775hpycHDw8POjTp8/fdnjR6dm6PFfqbKmTJk26yeFcHZMnT54UFxdH//79azoUwK5Cbdu2jfnz55Ofn09ycjJJSUmUl5dz++23U15eTkhICJ6ennTo0AEhBEIIgoOD+f7777Farfj5+eHt7Y3ZbGbBggX06NGDRx999KoVmTp16lC3bl1cHcRI6ijExsYSERHhTIXxK9zc3GjcuDGBgTdXEXB0/P39L1sHr8aYNw/W/YgidTRN50BqHsFnsqgoKYgqNVTstQVLhEZkbjkqEkXqGKREBRRdUCJKCc+zokhAShQpK7NWe1vLMZVbL5YhtHvAACnA3c8Mhfa8W0JKgsJ8CAnpi1/xBVbVqU2OVy2aFB7B67iBn4LaXPFUTooodokOHBONacMOAsllnejHLjoQoWYQSB5ubkVkF/gSXOxLpjmAqLICil1P42X1omVZbToqCYSlrKUgoB3PNmvCkrwCAkwGWtXy+KPegcsSGxtLZGQkwcHBdOnShf79+xMREcH06dOJjo7+Ww4vuri4EBkZ+ZsTyP6OzJ49m0mTJk2ubplDj7uYTCbOnz9fozGUlJQwffp0HnvsMc6fP89dd92Fv78/H3zwAR9++CGDBw9m/fr19O/fn7Zt29KwYUOSk5Mr/V/ffvsto0aN4j//+Q9JSUls27aNr776ildeeeWaze4XTflOLsWZQb4qzgzy1eOQGeQrClFLVcWmGlndqCO6yRVdVbGqBqyqEU1RUV1caTTqBWwuhgr7OUgh7DMPXVz5IeoRrAYFKUBVQaiAAGmAQ92CwSAq3PRgMQhsCpiNgqndhmEzGJFCQRoUsnonsyn4W8yl2xhcuIyDSnM+1cZSmH9tM663i1vx4gL15TFOi3CmiElMYwI6gvgWP5BfnESk5waE7k6d4jrkiiKaWO/gG/1xgijh0133EnYuj/a1PJiVdQabfvNHYaq7tzRu3Jj33nuPxYsXk52dfdNjqmmcGeSvD4dWtkaOHElGRgbx8fE1EsOUKVNYtGgR0dHRREZGEhMTQ3BwMBcuXOD7778nKiqKcePG4eXlRfPmzQF7eonAwEC+//57Tp48SZs2bejVqxeqqtKjRw9ycnLo3bs3jRs3vuZ4nMpW9TiVrao4la3qcUhlKyICAgIQ5eVYxz9Fxzf/hSkkCFFejnjmGWSfPhitZsTTT+P2xHjUwGAsJRcoDA3ngsWCGhzC9M4P8N/YnhhCQtAtGRzvZkQ2UDmLwkd9A5jStxZ53gbaFVkpvsXE1Hb++GuSzXfFMvnefxOmu9GowMKF9gPJ6WbFL/wgoWGHiXU7QGFxADt82rA36ur9rRc5IRqQL36pQXtahNOA4wSRzWkvyeGCIMweRdikRpsWG/AviMFaHou7xy4KrUai7nwef1cjX54+T4ynG408bu6973L3Fnd3d/r168fnn3+Ooij4V1fO6C+KU9m6PFdSthx60Ll+/fp8++233HfffTVy/JMnTxIZGcnatWspLi5m/PjxCCGIj48nOTmZzz77jJ9++omnn34aKSVbtmyhW7duDB48GLPZzJYtW+jRo8cl+xwyZMh1x/PJJ584PVvV4JyNWBXnbMTqccjZiNu3o49/CixmjJs2k1lYTtSrExFWDTUhASFAt1phy2Z7WoexYzFpGn6qipQSmZfH2PQZ9Bxbm9il0+xer1TFLmLpGs+m5qILwcT52ZisEiVd8IJShqpJYo/uZ59pMXfPfx+TVcM36wD1/adwvEcClnIjR8/XIcRmgau45XwkH2ZBzmOcVsNIDQzH02KlRU4eZw1WLrh5UGpyoczowtumF2goj3B35NfoGS4gof6ZpkQsXITWKYVGJUbqlifzXeBI6rkY6eVfi1AXIyvy8rkj8OZODrrSvcVkMvH+++8zevRoGjRo8LfxcDlnI14fDq1sjRgxgsTERAYMGFAjqkVqairPPvsscXFxdOvWjbKyskqDe2BgILm5uZw5cwY/Pz+WLl1K/fr1GT58OAAGg+GGz9ZwKlvV41S2quJUtqrHIZWtefOQ635E0XV0XedQajaR508jpETXNNA0FAmabuPI6Z8JyCpASAm6jtB1FOz+LHPWBWoXZNk9XFKiVni2VClodaEObudyUSUoUqDqEkWCLsGtoICoM7moUgcp0b0j8XAZjVdeHL6F3qyJjiKjlndluBPlqzSxJLHX0PaS0zBmuVAvpYD6Z0/TOuMIzbKOMTG7ET2yyhica4S8zTQ5lYafkkeadyRrRX/O+XlQq6gMowVcAvIprOWNX/E2aosTHEmNwXv7NrS8XPbW8ifVpvNweMBNfWt+694ihMDf358tW7Y43gzXPwinslU9UkrmzJnz5/RsAYSFhfHTTz/VyLHDw8MZMWIEZrMZT0/PKn6p7OxsunTpwh133MH48eMZN27cHxqP07NVPU7PVlWcnq3qcVTPllLh2VJcXQgY8QCayRVdESgqoApsQsFmMLCiuQvlJrvfyqaCpgjspabhrL8ZqahIAFVBU+w+LN0IWQPOIVUDulDQVRXdKNCEwKqYWB3ZGZtqxKqomI1GztazWyJcVXfCXevyZMalHVOR4k0HsZUH5aWqaWJoDN5WEz42N9xdC4mM3MnBBpMoilmHJWwbHh7nMUid6IwcpjKOB+R/OePuz9JWt7GidTDb3W4lsySTaMOPLLH1p1GBgeKEBHJffY2AFd+QVlzK+YRNN+EN+YWrube0b9+e9PT0mxOQA+D0bFXFZrMxd+7cK67j0MrW8OHDiYiI4Msvv6ycEXIzqV+/Pi1btmTUqFHs3bsXf3//yt68zWbjyJEjDB48GFdXV0JDQ//weJzKVvU4la2qOJWt6nFIZetXni3x9NP4Pz0OJTAQYTbDA72h1nGwanzXrgff1LuPAp9SvKSNxO6NKGvkTu3SEkpCvDDbXHEJtmHysJLZqgU/hw4hwOUMafcLcgZ4UObhivVsBD/fVZsFwQPAIvhv+54saHYnZ928qaPb+HLgA7w2sCd9s612dQyF2lZoHDaJH4x2S8RWv/ZstvRgp7Fj5SnUkvnkKiGUB52nq9cX1I05iJf3eczu5eQGFvB1YH0OhUaR7N+ACPcUGngdpiHH6C7WUDulBz+F+JAWFEx81m4CbGUUaf2pU7sXJbF1aPjeq5z39OIHowddNq+nTtfON+2tuZp7ixCCvLw8CgsL/xZqj1PZ+gWbzcaQIUMoKyvjhRdeYPr06ZdVthw6z1ZCQgJgT7swf/58Ro8eTatWrWokntGjR/OPf/yj8u8dO3bQtm1bbrvttpsWgzPPVvU482xVxZlnq3ocNc+W7F6RZ8vFBTFtGjw5Dqw2MBjsSpXVYp9d2McF1phBA3lxXEIHqSiUtXPDZVc5QtPRhYIiFBRpRTcITsYGEXk4F2GVWIUBEPb8XKqByd1H8sr6TzFqNqyqgfvv/TcvNtiNu/U2ZFkIXkZfjvYaTinuvMx/yBZh1Z6GiyzHLFypL48xJHEOGy6EsDfqQbRIb4xIAsUpSvCgDDf+xatkUodN3EY/viMxvxM7arXika0rCXGz0LDFCpRyE8UlrVhY927Wlf3yoL22VUOaet+cNBBXe28pLy9n/PjxNeYvvpk482z9gpSSMWPGsHLlSurWrYsQ4rJ5thx6GLFr164MGzaMTZs20aFDBz744AN27dp10+M4fvx4lV58SkrKJdngbwaPPfbYdc1i/KvjrI1YFWdtxOrp3r07t99+e02HcSkJCUizxZ6otNxM1mfzkBYL6DpYrfaOlgRdF5jzGiFt2McNNfuPkIAORceDUGwaqrQnOlU0K0IHxSbxyS6prH9o0DSMmhWD1DFqNgad+B4X3VL5d5eTe2kp1xCtTmSXdRZL0t/hQmID3CnlbZ7kHjnvkvC7HtlHdM5JzMKVrvJHThPGtCYvElLugyHcBIrgkeOljORDxmIf2p4kpvC5GE2KaMQH4hm2+LbHqFi4tfM8lDp7KdnzNq7HO/NmQF82l9amtdxZebyFh7ej67ab8tZc7b3F1dUVf39/ioqKqiwzm83s3LnT8WpyXifO2oi/sGLFCr744gvq1q37m+s6dGdr48aNPPvss/zwww+A/Qvkm2++uelxzJkzh06dOvHZZ58BkJeXR8OGDW/6sJXTs1U9Ts9WVZyerepxdM+WZjAyw7cpumrPlyWNRjC5oKsKZgPMjAvAbLCvK0wmMBnRBNhUlUMxd6IrCpqw98V0xYCuCHSj4PR9LkijsOfSUlWsqhGbULCqBr6p1xdpFOgKCJOkc9cdAPwgO/Bjswf57o7H+DzgCebzIAl0px3beb/gFW45kUpo/hlCCs/R7cQuvLULJIgePMwsgozZrOk1mIHr5lO7tIiZjTyZJP7D6+I1zMKt2sswmg8AqBeYymmP/7Jb8eO4Wp8RBxMYs3Mzc9Ke57PyYdxW8hgbfmjG7g0fU1hY+Ie+Nddybxk5ciQXR2N+jRCCLVu2sGTJElasWFFth+zPhNOzZWffvn3ExsbSrl27q1r/TzGMePbsWb777jumT59eI7Wp7r77bh5//HF69OhBQkIC8+fP59///vdNz63irI1YPc7aiFVx1kasHkeujcjSpTB4MNu6D0Kf+A86JW7k267j6deqIeq3y5hT38B77TMZtE7joeMWcm9ryk+nt9NpXxFnY0K5EPw8oXsWEZW+G1UYOBvQioBzp1DaZXBkrA8hKy8Qtrgcc8MHmGFS6X50G2sat2FB0ztYXvYo4dsKyL3VE5d4AxHp5Qxo9BUHQkKpVVaCRTVgNSrYhF3lEVKnUWEatc6Xo0hJfj2FI+KXWokBMpczIog35NOsYcAlRa1/TbQ8TB9WMZ2nieMgo5mOD/mc3jyMHYauLL3Fmzq5VvyKbYSVnuXWswdoUbQPH5GCSTtPlgwlO/4B+jz8xFWVPLtWrvXeMmbMGO6+++4q63/zzTdMmzaNjIwMZs2axdmzZ+nYsSNRUVE3POY/GmdtRPvkuN27d/POO+9c8vqVhhF/s7MlhHAFNgMu2PNyLZFSviKEqAt8DfgBe4EHpZQWIYQnMB/wBB6QUp4WQgwHvgCaSykPVOw3CRggpUy/zHHlr58SCgoKWLp0KdOmTcPX1/fKV+IGk5iYyOuvv05CQgJ33HEHgwcPZsCAATc1BnB6ti6H07NVFadnq3oc27NlQbiYENOmIZ8ci7BYL/VlqYK1g8K57ZssVJtEU0EgMGgSqaok9B9A55WrUTULulARQiJ0HYTg+Eg/Gsw5h7BIkKCpAkUHi2rgkzvuZNyqpQirRBoFR14N5PmId9kZHMvtluX0Sv+W84drcyHXnfyOjfiiyUNXfWrN5R4O0ozW7OIRPuFlbSo5hstP2nDXSxkuPiFyz0mOnh3EnkZ1OO/tQYGnByUVE4OC8200TTfTPKOQQPUgru6rOSeb0u6O+4ls2PCGfj9c671l3bp17N+/nzZtLi1rtGjRIj788MPKv8vKypg9ezZ79uzhoYeu/no6Ak7Plj156cyZM6t8Ln5vZ0sAHlLKYiGEEdgKjAeeBpZJKb8WQswE9kspPxZCjAYygZPA/VLK5yo6W68C26WUQyv2e02dLYCsrCzmzZvH6tWrrxjzH4GUkvHjxzNlyhQ8PG5+jS5wKluXw6lsVcWpbFWPQypbU6agv/hiZXHp0k6d8dySgJD22oVIe5UdXUBJQ1c8j5fbl1Vsbl+mkBnagojT+1CkrCznowA6UFDbm9oXLqDoVKaKUACbAmdauBGUWIaig67AD3fdxogx79JTX83j5z6j2EtFddEoyXVld0ZrPm33f3jJQorEHzejs8H5DLru/Qk/l4OE+ufi3aCEUoM7W2zd2So6k2WMREidjilJxJ9OQzelceu5DGpbrZz39KPbSwsxurn/7jiu9d4ipWTkyJE8+OCDl7y+dOlS3njjDTw9PS95/WIG+j+TSuRUtuyd57fffht390s/Y7/LIC/tFFf8aaz4kUA34OLA7RxgYMXvKvb/bx349UFXAnFCiEZXdzq/cO7cORYuXMiRI0d+M5fFH4UQgg8++KDGOlrg9GxdDqdnqypOz1b1OKpnS7i4oKsqFsXIv93iMKsmdIH9Dq3aO1qaQfBlW2/MBlHh07Ln2dIE6CpYW8aiq0Z0IdCFQFMr6h+aBAejvJAIdMCiUunRshoEh9p7VP4tjYKNcT1RpEZv80rqHTTj/4OgzpFyjG4aXdtupbdcRZGoRS2ZD8A/5RuVp9JF2/i7L0c96ylO+ISztntHom7NoXbDEgrywyk8GcWtJdt40/B/DDq1CCkUSgwCt+JwhLkem2s3wc/dk7qazvlnW6GZLb87lmu9twghaNq0KVlZWZe8HhYWRlJSUpX1hw4dyt69e393nDeTP4tnKz09nWXLlgH2TnBeXt4N23eTJk2uWfS5Ks+WEELFPlTYAPgQeBvYIaVsULE8AlgtpWwihKgNLABcsQ8tZlUoW62BXUB3KeWwq1G2Fi5cyIYNGwgKCuLhhx+mTp0613RyfzWcylb1OJWtqjiVrepxSGULYPt2SEjgbOsObPFvgLpjB7V2bSOrTggNzydQ+/hpSmKjEC1icD+dh5J6ikMdYjlVdpZ2R8toNmQ86i2dYPt2rD+sJ6VBC15Yl0jnkg/Z09id/Q3cGZ4dReD2XcT0sJeV8Ugs440wH0IauHNnRjHuWzUyAvx47v4PCbSeY5LXU2TuD+Hcjtoo6OhC0HzUEawYeJeJHBAt7fvRiylR7IpNkMzmDEF29UwohOWXcMrnyg+oL8qXeFu+hFkxVb7WWW5gB7fiV17CgK3fYapVRJ3QLAIDT6CqNt7Me4nk4Cbsbh/HwgVfUXQyE7fSMDyLflFbWseepckjrfHwuH4F5nruLUVFRUycOJGhQ4dWvnax2kh1w9fjxo1j0KBBf5r7159B2ZJSsmLFCjIzM6lbty7p6emcOHGCTp068cUXXxAXF8fo0aMJCAiook5dDbqus3z5cqZNm3bJ67879YOUUpNSNgfCgbbAZSuSSikLpJR9pZS3SSmz/mfxfKB9hd/rN5kyZQpFRUXk5uby+uuv869//YvTp08jpaSkpORv186cOZN9+/bVeByO1k6YMIGCgoIaj8OR2sOHD/Puu+/WeByO1n7//fesWLGixuOo0jZtinzuOdw6tmFg8zC6P9SXzl+8w8AJT9Dm/ZWELdlK86nLaHjXizSc8Bkh7y7l9mHv8NBDH9PivUWUN29RuR/DKy+QYvDG4BLMf5u9Q2KUfSLP7JB03hoUyOsBtVgQ5sWygf4cruNByE/h7EyOJcGvCTvr3U6hmztNjudw7nAtwuOzcQ8sQ0cBKTiT7INe6M2H896i1357ZY+LHS2AXBGCLhSksH+1XOxoGaX5svf593nmko4WwGbRjeF8QrZbbbzanKZjm/WEBR9FzQ3izJYHSPGLoZ2wkbV/H8VlJ0FIjLYivOjxx5gAACAASURBVERO5T5K8uaTdmLW73pfnn76aSwWyzVt5+npidFopLS0lLKyMqSUeHl5kZaWVu36vXr1Ys+ePUgpK9d35DYnJ4evv/66xuM4ePAgL7/8Mh9//DGLFy9m+fLlLFu2jJUrV7J8+XK8vb1ZsWIFQ4YM4Z133mHq1Kk888wzJCYm8uKLLyKlZObMmWiaVrnf8+fPk5iYWPn+XO74ZrOZsrIyCgsLL3k/r8Q1daWllAXARuxlSWsLIS5W3gwHTl3F9jbgHWDi1RzvzTff5LvvvmPQoEEsXLiQZs2aUa9ePWbPnk1ERARms/lv1Y4YMYJBgwbVeByO1r7++uvUr1+/xuNwpLZ3796MGTOmxuNwtHbChAn06tWrxuP4I9u8mTM4uyqTLko9Rp0KZsyOfxCU9MtMwTSzyroiIwvzXWh00gtZqILwILKlkYONmqPoOnWSt3DhpBdCAffAssptT20NJm31yxw73Z/WR7MBGJCWB1f4ohkmP+V9Hq92WVO5jyJRCzdZQrA8fcmyVdyJkDphZafITRzI8e/eJW37s2xz70i5yUSTjWc48mUWpkxf0BVMZfFomIjQDtA4ZCrytoNMmPDf33U9Dx48iK7r17zdnDlz2LhxI/fccw8Wi4WhQ4dy9uxZwsPDq6x/22238dJLL2GxWCrXd+TWx8eHrVu33vD9Tps27ZrWj46OZu/evQQHB/PZZ5/h4+PDl19+yWuvvcacOXN47rnn8PDwID4+nmbNmtGuXTvi4+MJDQ3lgQceYOjQoaxfv56kpCSGDBnCokWLuO+++wgNDWXQoEHs2bPnisePioqq8r5fiasxyAcAVillgRDCDVgLvAkMA5bKXwzyB6SUH11mH8OB1lLKsUIIE3AI8ALayWswyGdnZ7No0SKeeeYZOnbsWN1mf2mcsxGrxzkbsSrO2YjV45CzEW8kuckc6TWIEy17k+YysMpij7KTHK+9mu+aHUKq4GI28MDWfmjmAyBcWX3Ho6T61WbcvH/T+O6j6DaFo0vrIrVfnsuFIQwXr6Fsi8tjfZMYJhw4z07Xk2yObl5tSHt+KEIzlLA8qJivolzJcQvAJlQ0xVDt+heJ1E4yQFuKcZE7JhHFiWAfUqNCOBwRjFe5hdfWbMQgjewKKkNInX75p4gZ8xQl3kc4cvRFPNzr07r1ElT1+o3yv+feMnbsWDp16kRISAhgH36bP38+06dPr2IFmTBhAj179vxT3MNu1GzEHTt2UFBQQIcOHVAUhdtvv50xY8YQEBCAn58fwcHB+Pn5sXfvXpKTk3F1dWXIkCGXpPiw2WwYDAYSEhIwmUw899xz15T/sqSkhP79+9OxY0fGjRtX+V4BvPrqqwQHB9OoUfU2c4vFwrp163jrrbcqX7vSMOKVP+12QoA5Fb4tBVgkpVwphDgEfC2EeB1IBK7qri7t6SE+AN6/mvUvUlhYyOzZswH+tjXwHnvssUs+DE7sODPIV8WZQb56unfv/pfJ5F0tAY0piA2mb5056PocMvLj2VY0HHOEFYPrBcyFoYSXP8Bju91AKgQY1nDGvTs+ZnfyDbuIPXGYpJCunLkzlqaeyaSsbIHUyhFq8P+3d57hUVVrG75X2qRBSELoHULvoIgCgkg9AlIVlSOKh6gURRAFUTiKAoKKIAFBUYoCQvjoHHpROkiNgJQQSiCEFEhPZmZ/P2bAwKz0SWYyWfd15dozu76zWTx7vc9eBc1gekWn6W+Qln6cIzVbITSN6+5p1LkdyfHKicR7WLbP6tjeg7ZRbjhppdCn38PZ1UhqDioVMeml+Fn3H1JfdsdgfoninpZO7YgbNP/zMJdiTlHFuwEIHamu6ZQqV42LiVOIiziEn19bGjWcna+KFuRPW77++mtGjRpF8+bNqVWrFp6enrzwwguMGDGC4ODgh8aM7N+/P1u2bKFdu8Kb9zGvWGsE+WvXrgHQq1cvJkyYwOXLlylbtixXr14lLCyMS5cucfToUZo1a8aYMWM4ffo0M2bM4OWXX0an0xEfH8+CBQsoU6YMe/fuZdasWbmuG3h5ebFjxw7pGG0ff/wx77zzDiVLlpQ+d93c3B68RszJdYvEoKYAsbGxbNmyhTt37jBjxoxi2VheOVtylLNliXK25Di8s2Xm6sn9HN8ynV5JO0lzFfze2t9in4hDr3Mv/AnudxrXNA2D6y8se74HV5yrUzYtlmGrz5OS7oyzW00ADGmX0KccxGiI5FLDNhxu9AQ3fErhZDRgdMrZoKJV79yk7L0YKsbd4aaPP4dqNKA816kbdYX6AUc4dfdxwtMr4Ww04O96jzI+V/EQSZS+5I3fidaUrhpN+DUNt6il+FdqzTmfeNJFGgMTn+FX712U1JegWquWDO6W/7EQ86stRqORjz/+mICAAJo1awaYxoxcvXo18+bNw8PDNJq+pmkEBQXx8ssv5zvmgsZaztaaNWsYPHgwDRo0wMUlJ76PafincePG4e/vj9Fo5Pfff+fYsWMMGTKEmjVrEhERQXCw9AVbnkhLSyMoKIg+ffpQsmRJi+179uyhV69eNGnSBMjnOFu24tHKll6v59tvv2XGjBnUrl3bdoHZENUbUY7qjWiJ6o0ox257IxYAxy8dZMdHn2Pw8ManfjqVq4fj6x5Fsuc/laKUq7W4diIIQ0opdKWuUq3dNFJw4S3dQtKdXJm08jcMxmcfOq+maTT+cxzhgfXpVEnPPp9BrK6osaNquYf2+zTiJ1pHbmSPXytW+3Tl7xL10DubHqpCM6JLTyfFTUfThFDe8ZqMlijQ48KqS13xTINySSZXytU1mWq1f6ecfyRJUb5cPjqQq+IcaZ4GhJsvfmmWs4r8iQfLxo3GS5ezh3hmWEtbvvnmGxITE2nbti1gmhVl06ZNBAcHP6jITZ06lYoVK9p92cxvb8SYmBg2btxIv3796Nq1a66Pv3fvHlevXqVhw4YkJibi5eXFzp07+eqrr/j8889p2lT+Ojuv3L17l2HDhvHaa69ZVAoTEhI4cuQIkyZNAhyksrV//350Oh2jRo0qtq8RlbMlRzlblihnS05xcbbuc+Lmnxw7upOIEydx+TsOo28lRHkvHmu5DidnA04GDQxO3Al7Fp9q+zDqdQQeSmRB+c78XP812ly5zsgf57Ct+zuUjHfCM/afh42LPpGqaceoV7YtmnsK93wuMKpRacJdq/I9gy1iScGdP1OfIjK1Gnc0bxKdvAhIjWNg6bk48c9zKDa+Mul6D6JiPLkbVQ1DWglAo3TAFWoFHiIxwZfTp7qgF3qSnVMpof/n1WWc6z2O+XdmYrtmPFu/bL7vnzW1ZenSpZw8efLB7CO3bt1i586dzJkzBxcXF/R6PcOGDePZZ5+lTJnMR9m3Nbl1tq5evcqmTZsoUaIEt27don79+owbNw4fH+sOipvT13l54eOPP6ZNmza4m2cxyMjy5cuZO3cu4CCVrcjISM6fP8/HH39su6BsjHK25ChnyxLlbMkpTs7WoyQnxLN21VwuHL/IY93/4E5UFfrcOMSh5r4YnU3Ph0sbPyftXhhrAnTQoCpnKtTE2WjA4OSM0DSGHYnFN8yIV9m/cPe/grtvOJ6+4Th7xrKdLvwkhtIt9i+G7ovjrqsH6V53SXGLI83tHpp3JJ4+EXj7ROLkkorsuZhyvRF43sXVOwpnN1MvyJR0N+KSPLh4rRH+3jHUq/YX+7Y1o2TMy7gKH8JLHyZFF0PVyKcwuvjwZ213pg1sQhVf+YTXucHa2rJx40YOHDhAp06dANNrsdDQUL74wjQwbEpKCkFBQfTr10/62soeyK2zFR0dzYwZMyhfvjxDhgyhe/fuBRyh9XnnnXfo06ePdNumTZt4++23qVq1atGtbG3bto3Dhw8TERGB0WjExcWFWbNm2To0m6GcLTnK2bJEOVtyipuzJUOfns6uPU+SlGQkIbwKFeucwOBiqkzEX26Od+J1RKPbpOlh+fn3iPSpSLWIq/zeoAkxXiVpm3iCF71mUVKLJy2+LCmxVfjDpTG/VuxK4zsJfLzvGiWFM+4uJXF3cideH8vVhL+IS7uNAVechSdp6Te4a7iJ5ipILV+XF107E53szp+aE43aV6JOO39Wbf0aw5W9VCkfiWu5NBCCe1dbUrLqYfTJLpT/riZVpn3D1u9OcsfZNJZYgqcT3klG4mt58eGYVvm+V9bWlpCQECIjI2nQ4J+hOFasWPFQO6O4uDiGDx/OoEGDpE6KrcmpsxUXF8f69eupUaMGr732WpHu3PXWW2/x4osvSrdFRUURHh7O6NGji25la8uWLXz77bfMmTOnSM6Obm2UsyVHOVuWKGdLTnF2tjISFbWNU6ffBowYkl154vQdTlcqT0qFJESMwP13Z1zD3DB4GPA5qXGteQ0O+rTgWNtAjparR0lDIiM37iItrgTCuTRrO5bjbNkAPloTh/7eAfQpB9CVGgVaCjUvLET4leeuvxth+kSSDUl4GePw9Esk/dmS/F5uKgN37aJF526c2hXLjfOxlK7sTYM2FfD2M1U2dv96kFL1fqREpVMAGFKcqPyeC7p69cDrKSLLNee4wQlDioFkN4GrXuOZzx6niZ93Fnche6ytLUOHDmXgwIEPnW/NmjV88skn+Pr6PlgXERHBmDFjKFu2LJUrV6ZFixZWub41yKmzdX+csfr16xdSZAVDZGQkL7zwAo0bN6Zv377Sfe5PNJ7vEeRthU6no3Llyg8VwuKMmhtRjpob0RI1N6Icu5wb0QYEBHTiqSf30LDR9/zt1I9N15vz2MUI6lxIQPM1kt7JiGeXVFrUv4WTzsi9BCdW9OhAf/flfKL/mBSjjmndu5NQXk96wkpa3byD0cmJc5670KeeMF/FiHDyoIIxmlrRW7n51CXmlenJ/2o+Q4vX/qJWz+ukVh/Blw2q0X/sW9RoWo1e7zaly38akpKYzp5lf7Nxzik2zjmF0HypeHIEV49VITHKk4S1JpckLTwcl9KBBDYuzWuftabO4+XwStVwM8Dti3fzfZ+soS2hoaF06dKF0aNHEx4ezrx58x7a3qBBA4s5fytUqMCvv/7KpEmTOHfuXL6ub21yOjfirVu3qFs308lmigx+fn5MnDgxyx6Trq6uxMTEZHke5/ut6O2N//73v5M6d+6Mk5MTHTt2tHU4dkHVqlWpXr26XVrLtqR+/fpUrly52HackOHh4UG9evXsuqGtLShdujSBgYFWb5xbFHFxKYG3Vw0eD3yG6/7+hFyM4HZqUzSjOzqfGBLKGjnpV4PzPjWo4/YX7WttoIR3FHE3y+N13sDVMpU4EFibSs4BPHnsINsaNyLVx58mZy+Aex183XVoBkGUR11KJoXzhfdgUo1ONPTdSlhCGp3bL6R39Sfxc/3nISaEwK+CF407VKZBm4rUalGG6o1L82TfWrhfTaCW4Rk8v1pNQHgiAB6NGuFa9SmEsxulnqpIzWYBlK/lw707yTz2WHm8S+Xv9Z81tKVMmTJ4enpy4cIF/vOf/9CyZcuHnC1/f3/Onj3L7du3H3q9CODu7s7u3bszHVjTFuh0OqpUqZKlCRIfH8+lS5cedAYoyjg7OzN//nyefvrpTF8ne3p68tdff7Fu3TomTZr0X9k+du1s7dy5kxEjRtg6DLtBOVtylLNliXK25ChnS87/lm3l6c4f0eSWK973rpJicOXmtUqU8Iqk+tMnuNfHgHu4kVK/OFPpj3jq3LpBp1PLKJ14j+UtmjDvX91ofe0Gf5erwPVaT+Dq2YGAa8eocnkLCSWq8F7DNzGIdEZWP8rkY2l88HsN6gU0yzQeJyeBt6+OcjV8qN4kADcPFzyblUF/OwX3lm+AkwulXniBgPe+xpjohJP3P3MrVq7rR58xLShbLf8NzK2lLX379iUoKIiffvqJ5ORki+0dOnTg8OHD7Nixw2JbYGAgO3fuzHcM1iI7Zys+Pp6XX345W6enKPH666/z6Iw2GalSpQrHjx/P8hx23WZr9uzZDB8+3Nah2A2qzZYc1WbLEtVmS45qsyXnvrZ4xMZy8qcexLdNplXT9bj7ViYmdj+eHtX5+9OX0G1JwCnV5PIYBNzz9WbBy6PYUrcuOn0ayW7u1Lhzh9G/HOZWuSeobvyDMKc2RLv/zWPtvena62UixrxP8qlT1Nq2NVcxGpKSCB/8GW41uiM8BLpqpUg5F4tbZU9KdiyJLrAmwsoaYG1tuXPnDmPGjKFTp05cuXKFs2fP8u9//xudToemaaxYsYI33njDoo3Whg0bWLduHQMGDLD5bBk5bbMVEhLiUB3aPvroIyIjI2nSpAmNGze22L5ixQrmzp1bNNtsDR061NYh2BXK2ZKjnC1LlLMlRzlbcu5ri656dWp1nYZIgRM7XyD13CXKBHQhNeoasV1jMQ6vz4FPRvNVbyc2tgxAV6Ya/z35ByN2raNMkmki3gCDO/XOLcXr7hHCnNoA4J9Sm/ZPPoUQAtcK5Um/dQtDQkKuYkw6dIjUU2tI/Wsx+htnSD5xlfTrx4j+bjBhvXpyd906q98Xa2tL6dKlWbBgAefOnaNZs2ZMmzaNpUuXotfrEUIwYMAA5syZw8WLFx867rnnnuODDz5g0aJFREVFWS2evJDTNltpaWnYq5mTFz7//HPmzZvHX3/9Jd2eXfMeu3a27DU2W6GcLTnK2bJEOVtylLMl51FtCd8yjYuu83E/6kSFhE7canKCFOdInmi2Ga18VZ5Y/CxG5wSc9SVxNrqT5haFX1ogj6UPoOeu3ZQJ3Q7A7YC6nGlgagriJFJ4slcatSvW4kr//vi/9SZl3nknxzFq6enErVpF4v79CHcPXPx8cfb1Je36de6uCqHKzz/h9cQTVr0vhaEt586dY/r06QwaNAghBHq9nlWrVj0YKDMjSUlJfPDBB9SpU4dGjRoVWExZkVNna9u2bQwZMoTq1asXUmSFw8iRI6W9Eou0s6V4GOVsyVHOliXK2ZKjnC05j2pL1S4fUK3Cm6S00LjcfgtJvpGUDquPR9XaeLrpmNhkCjW09viJWuicfKnO03So8Rx9ezTiuttd0p2dca5UiTJR52h8ZSIARs2dsD1L8WhQn5LduxPz08+kXr4sjSf99m3id+8mKjiY6yNGcLnX89yc8DHC3YOy48ZRcfqXlB03jtJvvokxIREnb288zPPTWZPC0Ja6desSFBTE6tWrAXBxccHNzY24uDiLfT09PZk5cyZHjhwp0JiyIqfOVrVq1Th69GghRFS4eHt7k5aW9tA6vV6Pl5flJOwZUc5WEUI5W3KUs2WJcrbkKGdLTmbakpJ6i/O/j+Xulf00b7oU75bZO0dJUVFc+O9nuO3cwd1Br3LB+Tx1o9IIT69PmYrnaTdqPmk3I7n8XA+01FTcatTA66mncPL0JOXsX6T8dRbDnTsPzudWrRqulSqREhqKITYWANeqVfB6vBUu5cpyZ/Z3+A8dSpn3Rln3plC42jJ69Gh69uwJmAYOFUIwePBgi/1WrFhBdHS0Rc/FwiKnzlbGwT4diW3btrF+/Xp69er1YDiI0NBQKlSoQO/evYvmoKb2GputUCPIy1EjyFuiRpCXo0aQl1MQ2mJMSsLJ0xOS4wjbv45KLbvg6vPPfIVpV68Sv3Mnifv2k3TkCJpej65mTdzr18e9Xj3cG9RHV6cuzt4mx0AzGkm9cIGkQ4dIPHiIpMOHMWZo9+XRrBkIgVv1apSfPNkqQ8EUpraMGTOGHj16AKZ5/kJCQpg9e7bFfu+99x4VKlSgZcuWNhnuJqcjyG/cuJFhw4ZRtWrVQoqs8Dh+/Dhz5syhUaNGNG3alJCQECZPnoyPj49jVLbCwsIc7v1vblDOlhzlbFminC05ytmSY2ttMaalgabhlItKjabXk/LXX1wZ8ILFtrqnTyGs0GuvMLVl7Nix/Otf/3rwfcWKFXz11Vd4eDw8x6OmaWzevJmQkBDq169PixYtCrXSlZ2zlZqayvbt21m6dCkRERGFFpctePPNNxk4cGDRH0E+Izt37mTcuHG89957pKam2jocm6DabMlRbbYsUW225Kg2W3JsrS1Obm65qmgBCBcXPBo3pua2rTh5eyPc3dHVr0fZcR9apaIFhastLi4uGI3GB98bNWrEpk2bLPYTQtC9e3d++OEH6tWrx5IlSyx6LxYk2bXZ2r9/PzExMTz++OOFFpMt0DSNixcvsmDBgkx7KGakSDhb165dY+LEifz73/9m3bp1jBw5sljOlWjr7NNeUc6WJcrZkqOcLTlFXVuMqakIV1e7H2crK6ZMmUJgYCABAQEAGI1G1q5dy8yZM7M8TtM03n777UwnSrY22Tlb69ev54033qBSpUp4e+dvbsqiQHh4OGFhYbRv375oO1upqal89NFHvPjiiwghSE5OLpYVLbB99mmvKGfLEuVsyVHOlpyiri1OOp3VK1pQuNry9ttvs337dkJDQwFwcnIiJSUl2+sLIUhJSZGOTF8QZOdsJSYmEhgYWCwqWmCaRq99+/bZ7mf3la3x48fTvXt3dDodUVFR1KxZ09Yh2YygoCDq1atn6zDsjunTp9t8VGV7o2rVqowaZf3eWUWdjh07PmiErPgHpS1yClNbfHx8HoyttXbtWoxGI40aNWLs2LFczmSIjPtMmTKFX375hYRcDhSbF8qVK0e/fv2y3MfZ2bnA4yhq2HVl64cffqBSpUqUL2+a4f3AgQO89NJLNo7KdhT17LOgUM6WJcrZkqOcLTlKW+QUtrYIIXj77bcZNGgQP/30ExUqVKBjx458//33BAUF8ffff0uPK1euHDNnzmT58uXcvXu3QGPM6Thbioex68rWX3/9RfPmzR98T0pKokKFCjaMyLao7FOOcrYsUc6WHOVsyVHaIsdW2tKsWTPmzJnD9u3buXTpEt27d6dv377Mnz8/02P8/f2ZM2cOISEhREdHF1hs2TlbthiOoihg15Wt+wO8galha3Fv7KuyTznK2bJEOVtylLMlR2mLHFtqS4kSJfjuu+/w9PRk165d6HQ6EhISsmybVbJkSebOncvGjRu5fft2gcSlnK28YdeVrYw15MOHDxdabwt7RWWfcpSzZYlytuQoZ0uO0hY5ttYWIQT/+c9/iImJQa/X07JlS5YvX57lMZ6engQHB7Nt2zauX79u9Zhy0mbL0UlPT+fUqVOcOXPmoeE6ssKuK1sZSU9Px9/f39Zh2BSVfcpRzpYlytmSo5wtOUpb5NiLtgwaNIglS5awb98+WrVqle3+Op2OOXPmsG/fPsLDw60aS3bOlqO9RoyOjmb79u1MnTqV9957j+HDh/POO++wZs0aQkJCGDZsGO+++y7Tpk3L8jx2Pc7W/PnzqV279oOeGK+99pqNo7ItRX0snIJCjbNliRpnS44aZ0uO0hY59qQtW7du5dlnn81VLAaDgTfeeINXXnnlwTx++SW7cbZWrFhBcHCwVa5lCyIiIpg9ezbJycmkpaXh6urKjRs3qF69OqVLl8bHx0fqAsfFxfH8888XzXG2Dhw4QHJyMseOHSv2tiWo7DMz7CX7tCeUsyVHOVtylLbIsSdt6dy5c64rfc7OzowcOZJdu3ZZLY7snC2j0Yi9mjg5Yf78+Tz55JP07t2bJ554gtjYWN59911effVVKlSo8MBhfJTsEhW7rmxNnjyZVatWIYSgRIkStg7H5qh2FXJs3a7CHlFttuSoNltylLbIcQRtadasGdHR0VarNGbXZisgIKBQpw+yNrdu3cLLy4v169dz48YNFi5cSJs2bWjYsCE6nY5y5crh4+NDUlJSrs5r15WtChUq8PTTT9O1a1dbh2IXqOxTjj1ln/aCcrbkKGdLjtIWOY6iLW+//TY7duywyrmyc7YaNWpESEiIVa5VWJw/f54BAwZQv359Ll++zKJFixg0aBDjxo176PXrgAED+OGHH6hWrRqHDx9mzZo1/Pbbb6xevZojR45keQ27brNlr7HZCtWuQo49tauwF1SbLTmqzZYcpS1yHElbhg8fTs+ePdHlcsLvR8muzRbAsmXLmDp1qs3KU1JSEvv376d69erUqFEj20b7Z8+e5ejRo9SrVw8vLy9q1aqVK0czPj6ekydP0rZt20zbbKnKVhHi3Xff5YUXXqB169a2DsWuaNeuHdu2bcu3iDgSZ86c4ZtvvuHHH3+0dSh2xeLFi4mPj2fYsGG2DsWuUNoix5G05Y8//mDXrl20a9cuX+e5fPkyq1atYuzYsZnuExMTw7Fjx5g8eXK+rpUbLl68yNq1a7ly5QqaplGzZk1iYmKIiorC3d2dsmXL0rdvX2rXrl1gMWQ1EXW2lS0hRGVgMVAW0ID5mqZ9K4TwA1YA1YArwABN02KFEE7Az0At4D+apoUKIdoDu4CemqatN593AzBD07TdmVxXVbYeQWWfchwp+7QWytmSo5wtOUpb5DiStowbN46nnnoq3+2fc+JsAaxcuZJx48ZRsWLFfF0vMzRNY+vWrezYsYP4+Hj8/Pxo3rw5pUuXlu5/9OhRfH19ef311wskHsi6spWTEqQHRmuaVh94AhgmhKgPfAjs0DQtENhh/g7QGTgE9AZGZzjPdeCjvP0EBah2FZnhKO0qrIlqsyVHtdmSo7RFjqNoS0xMDPfu3bNKR7OcjiDfvXt3vvrqq3xfLzNOnDjBpk2b6NKlCy+++CKdO3fOtKIFEBoayquvvlpg8WRHtgNvaJp2E7hp/hwvhDgLVAR6Ae3Nuy0CdgMfAM6A0fyXsYZ3EnAVQnTSNG2bleIvVgQFBT2YlFvxD47QY8jaqN6Icjp27IjBYLB1GHaH0hY5jqItqampXL58mXv37lGyZMl8nSunI8h7eXnh5OTE+fPnqVOnTr6uKWP79u106NAhR+OH/fnnn/Ts2RNnZ2erx5FTcuWNCiGqAc0wOVdlzRUxgFuYXjMCbAGeBtYBXz9yis+BCXmMtdijsk85jpJ9H5Tz9QAAIABJREFUWhPlbMnJj7N148YNK0djPyhtkeMo2jJx4kTatWuHl5dXvs+Vm7kRO3fuzLfffpvva8q4ceMGvr6+2e5nNBo5ffo0vXr1KpA4ckqOK1tCCG8gBHhX07R7GbeZG1dp5s96TdNe1DSttaZppx/Zb6/5XG1yck29Xs/nn3+ulublkCFD2Lhxo83jsLfllClT+PLLL20ehz0tf/vtN0aMGGHzOOxtGRoaSrdu3fJ0/LBhw3j66adJTk62+e+w9lLTNAIDA20eh70t69evjxDC5nHkdxkREUGZMmX49ddf0ev1LFmyJM/L0qVLo9PpcrT/b7/9hp+fH2+++aZVf8/kyZNJSkrKUbx79uzB2dkZg8FQ4Pc535UtIYSruaL1i6Zpq82rI4UQ5c3bywM5nWI8x+6WpmkYDAa1NC8XLFhAZGSkzeOwt+WHH35Ienq6zeOwp2V0dDTffvutzeOwt+WFCxfYsGFDno6vWrUqAQEBDB06lNjYWLv4PdZaHj58mLNnz9o8Dntb/u9//yMtLc3mceR32aJFC44cOUJ0dDTAg8mT87KMjIzkzJkzOd6/Xbt2HDlyxKq/59atW5QtWzbb6+v1em7cuEHlypUL5T5nRU56IwpMbbJiNE17N8P66UC0pmlThRAfAn6apkn7gpp7I47RNO058/dDQHng35rqjZhjVI8hOY7UY8haqN6IcvLTG3H06NH07NmT1NRU5s6dy3//+18aN25cAFEWPkpb5DiSthiNRsaPH0/FihXzVW5z2hsxIwcPHqRJkyZWG6B89uzZlCtXjjJlymS5X1xcHCtWrOD111/n2Weftcq1syK/vRGfAgYBzwghTpj/ugNTgU5CiAvAs+bvOeVzoHIu9leg2lVkhqO0q7Amqs2WnPy02XJyciIuLo7Vq1fTpk0b6tata+XobIfSFjmOpC1OTk5MmTKF5ORkDhw4kOfz5KbN1n1atWrFihUrHjhO+SElJYUjR45kW9EC03yFNWvWpHJl21c31KCmRQiVfcpxpOzTWihnS05+nK33338fJycnxo4di7+/fwFEZzuUtshxVG358ccfCQsLo1OnTrk+Ni/OFpiGXihZsiSvvPJKrq95H6PRyLBhw3jmmWeyrWxpmoYQguXLlxMcHJztKPLWIL/OlsJOUNmnHEfKPq2Fcrbk5MfZmj59OtOmTXO4ihYobckMR9WWIUOG0KJFC/7v//4v27ZGj5IXZwugQYMGbN++nbFjxzJ27Fj27t2b5f5Go5F9+/Zx6NChBzF++umntGjRIkeu1i+//MKhQ4fw8vIqlIpWdti1s3X37l0GDx5M5cqVcXFxwcnJiYCAAKpUqULHjh0JCAiwdZiFiso+5Thq9pkflLMlR40gL0dpixxH15a9e/eyZMkSBg4cmOMxqPLqbME/bpOmaaxbt47WrVvTv3//B9vT0tLYtm0b27dvJzk5mWrVqmE0GgkPDyctLY369evz+OOPZ3ud0NBQLl68yMGDB/nwww/p2bNnrmPNC/marsdWCCG0v//+m19++YUOHToApn+o+Ph4duzYQZ8+ffI9x1NRQ81fJseR5i+zFmpuRDlqbkQ5SlvkFAdtCQ0N5ccffyQ1NRV3d3fq1KlDzZo1Mx0sNCdzI+aUXbt24efnR+nSpTl48CDp6enUq1ePBg0aWFz/fkUtJ9y+fZtFixZRtWpVpk+fjqenZ75jzQlFtrL1xx9/cODAAR577LGHtv3666/MmTMnRyPHOhIq+5Tj6NlnXlDOlhzlbMlR2iKnuGlLVFQUu3fv5siRI8TGxvLcc89ZlIn8OFsyQkNDcXNzo1atWlZ93bdixQqCg4Otdr6cUGTbbN2+fRsfHx+L9Zqm2XTYfVuh2lXIcdR2FflBtdmSo+ZGlKO0RU5x05aAgAD69+/Pl19+yddff82GDRsICwt7aJ+8ttnKjAYNGhAYGGj1dlX2VkG2r2ge4datW9LKVpUqVThy5IgNIrItQUFB1KtXz9Zh2B2OMn+ZNVFzI8rp2LEjPXr0sHUYdofSFjnFWVtKlCjBvHnzuHz5MgcOHHjQSD2ncyPaGlXZygWZOVuPP/44ISEhNojItqjsU05xyz5zgnK25ChnS47SFjnFXVucnZ359NNPadCgAStXrmTNmjWcO3fOqs5WQWFvzYzsK5pHSExMxN3d3WK9m5sbUVFRNojItgQFBVG+fHlbh2F3FOfsMzOUsyWnY8eOGAwGW4dhdyhtkaO0xUSPHj3o0aMHkZGRzJ8/Hw8PD/bu3UurVq3ssvNAeHg4FSpUsHUYD2HXzlZmo81u376d119/vZCjsT0q+5RT3LNPGcrZkqOcLTlKW+QobXmYsmXL0rt3b0qWLEm3bt3Ytm0by5Yt49y5c7ker6ugMBqNbN++nXfeecfWoTyEXfdGHDhwIEFBQQ+tv3nzJmfOnGHy5Mk2isx2qB5Dcopbj6GcoHojylG9EeUobZGjtMWSR7UlNTWVNWvWsGfPHnQ6HW3atKF06dI5Pt/9SZyt9dpv8+bN9O3b12IUg8KgyPZGfHTiSKPRyKZNm5gwYYKNIrItKvuUo7JPS5SzJUc5W3KUtshR2mLJo9qi0+l44YUXCA4O5oMPPuDixYv88ssv7N27l9TUVIvjk5KSOHbsGCtXrmTlypVs3ryZxYsXc/369XzHFhkZiRAiy4pWWloaV69eZf/+/axYsYIZM2YwadIkrl69mu/rZ4VdO1u7d+9+aN3mzZvp06dPjkaQdURU9ilHZZ+WKGdLjnK25ChtkaO0xZKcaIumaRw9epQVK1aQkJBAQEAAUVFR6HQ6/Pz8aNu2La1atcLLywuA9PR0Jk2ahKenJ23atMlzbLdv32bPnj0A+Pv7k5aWRnp6OgaD4cEfgK+vLz4+Pvj7++Pv74+zszM7d+7EYDAwdOhQGjRokKfrF9lBTTNWtq5fv87FixeZOHGi7YKyMWqUZznFYZTn3KJGkJejRpCXo7RFjtIWS3KrLampqVy5coXAwMBsK60hISFs3bqVfv364ebmlq84k5OTcXNzy9WYnKmpqezevZvY2FhefvnlXFf8HKKytW/fPp599tli62qByj4zQ2WflihnS45ytuQobZGjtMWSgtaWsLAwJk6cSLdu3WzWo9BgMLB//350Oh3vv/9+jo8rsm224uPjiY6OBqBp06Zs2rTJxhHZFtWuQo5qV2GJarMlR7XZkqO0RY7SFksKWluqV6/ODz/8wMmTJ9m/f3+BXScrnJ2dadu2rVXakd3Hrp2t/v37U6VKFdzd3enUqZNN5jqyJ1T2KUdln5YoZ0uOcrbkKG2Ro7TFksLUlpUrV7Jz50769u1rk/HOli9fzty5c3O8f5F1tpo3b07nzp3x9fUFoEyZMmzatInY2FiWLVvGlClT0Ov1No6y8FDZpxyVfVqinC05ytmSo7RFjtIWSwpTW/r378/777/PkiVLuHnzZqFcMyMeHh7cuXPHKueya2dL0zSGDx9Onz59cHZ2Ji0tjQ0bNnD27Fk6duyIp6cnHh4eXLlyhZYtW9K9e/d8N6qzZ1T2KUdln5YoZ0uOcrbkKG2Ro7TFEltoS2pqKhMnTqRkyZK0bt260P49li9fzuzZs3M8BliRdbaMRiMXLlzg9OnT3Llzh8WLF3P37l1atmxJQkICR44cYcuWLbi6unLjxg0GDBhASkqKrcMuMFT2KUdln5YoZ0uOcrbkKG2Ro7TFEltoi06nY+rUqbRs2ZL169ezbNky9uzZw7179wr0uq6urlYbbNXuna309HR27NjB/v37GTt2LK6urowYMYISJUrg6uqKwWCgS5cuzJ07l4EDB9KnTx+HzUJU9ilHZZ+WKGdLjnK25ChtkaO0xRJ70BZN0zh9+jQbNmzgxo0bODs707BhQ4vhJRISEoiMjMTX1xc/P79cXeN+3WPKlCk5PqbIOltgqll27dqVTz/9FG9vb3Q6HXPmzOHq1au0bt2aUqVKceXKFapUqVKo9qItUNmnHJV9WqKcLTnK2ZKjtEWO0hZL7EFbhBA0btyY8ePHM2fOHD777DN8fHxYt24db7zxBv3792f9+vWcPHkST09Ptm7dmutrXLp0yapDTdm9s5UZBoMBvV6Pq6srEyZM4M6dOxw9ehRN0zh+/HghRlp4qOxTjso+LbGH7NMeUc6WHKUtcpS2WGLv2qJpGjdv3qRChQpomsb06dPZv38/o0aNytV51q9fz/vvv0/ZsmVzfEyRdrYyw9nZGZ1Oh5OTE1988QXz589nzZo1VKtWjc2bNz/4u3v3rq1DtRoq+5Sjsk9L7CH7tEeUsyVHaYscpS2W2Lu2CCEeDIYqhMDFxYWWLVvm+jxJSUm5qmhlG1dRdbZkaJrG4cOHGTt2LL1798bZ2ZnDhw/z/fff4+npWUCRFh4q+5Sjsk9L7D37tBXK2ZKjtEWO0hZLiqK2/Pjjj4SHh9OxY8ccH7NmzZpcVyod0tmSIYSgVatWLFq0iBMnTtCwYUN69+6dqwZu9ozKPuWo7NMSe88+bYVytuQobZGjtMWSoqgtQ4YMoWHDhmzYsCFH+8fGxlp9qqAi52x9+umnlCpViqFDh+Lu7p7p8cePH+eHH35gwIABhISEMGvWrIIMt1BQ2acclX1aUhSzz8JAOVtylLbIUdpiSVHWls2bN7Np0yb69u2LEFIDCoADBw7wzDPP5LqBvEM5W/Hx8fj6+jJixAi+/fZbJk2axMmTJy32a9as2YMBTtPS0jAajYUdqtVR2acclX1aUhSzz8JAOVtylLbIUdpiSVHWlm7dujFgwAB+/fXXLOsEV69epWnTpla9draVLSHEQiHEbSHEmQzr/IQQ24QQF8xLX/N6JyHEYiHEfiFEA/O69kIITQjRI8PxG4QQ7fMSsKurKxUrVuSVV17B39+f6OjoTAcy1ev17NmzB19fXy5cuJCXy9kVQUFB1KtXz9Zh2B3Tp0+3ybxZ9kzVqlVz3fumONCxY0d69OiR/Y7FDKUtcpS2WFLUtaVt27a8/fbbLF68ONPp/lJTU602mOl9cuJs/Qx0fWTdh8AOTdMCgR3m7wCdgUNAb2B0hv2vAx/lK1Iz/v7+xMfHo2kaly9fRq/Xs379epKSkiz2nT17Nt27dycuLo5bt25Z4/I2RWWfclT2aUlRzj4LEuVsyVHaIkdpiyWOoC1NmzbFz8+PY8eOSbe3bt3a6r8x28qWpml7gZhHVvcCFpk/LwKeN392Bozmv4zvLU8Cd4UQnfIVLabJqMPCwvj555954oknCAwMJCAggDFjxkj3b9myJXPnzuXpp5/O76Vtjso+5ajs05Kinn0WFMrZkqO0RY7SFkscRVsaN25MZGSkdFtgYCARERGcOHHCatfLa5utspqm3Z+C+xZwfzCKLcDTwDrg60eO+RyYkMfrPaBatWocPnyY7777js6dO3P+/HmSk5MZMGBAfk9t96jsU47KPi1xhOyzIFDOlhylLXKUtljiKNry6quvkpiYmOn2f/3rX3z99ddWm2853w3kzV0GNfNnvaZpL2qa1lrTtNOP7LcXQAjRJqfnNhgMhISEPLS8ffs2S5YsYcuWLezbt4+oqCjOnz/PnTt3pPs70vKNN97gwoULNo/D3pZTp05l3bp1No/DnpYnTpxg5MiRNo/D3pYpKSl0797d5nHY27J69erUrl3b5nHY27Jbt244OTnZPA57WlaqVImGDRvaPA5rLC9dusTmzZsZP348MTEx7NmzB4PBwJ49ewAoXbo0EydOzPH5siKvla1IIUR5c+WpPHA7h8flyt3K6keNHz+eNWvWEBcXR4kSJVi9erV0v19++YVDhw7Z/B/VGsvvv/+epUuX2jwOe1t+8MEHrFq1yuZx2NPy119/ZebMmTaPw96WP/30E+vWrbN5HPa2DA4OJjQ01OZx2Nvym2++ISUlxeZx2NMyLCyM7777zuZx5Hep1+s5d+4cLVu25Pr168yfP5/169djNBrZu3cvRqORU6dO4ebmluPfmxU5GmdLCFEN2KBpWkPz9+lAtKZpU4UQHwJ+mqaNzeTY9sAYTdOeM38/BJQH/q1p2u4srpnlCPLBwcEYDAYaN26c5UivZ86cYeTIkbRv355x48YV6ffvaiwcOWosHEuK8lg4BYkaZ0uO0hY5SlsscSRtiYuLY/LkycTExDBv3jxmzZpFYmIi7du3f7CPpmksXryYr776Cl9f3yzPl69xtoQQy4ADQB0hxHUhxBBgKtBJCHEBeNb8Pad8DlTOxf5S3nrrLS5cuMDNmzeJi4vLdL8qVarQpEkTatSowdChQ/N7WZui2lXIUe0qLHGUdhXWRrXZkqO0RY7SFkscSVtKlSrFjBkzWLhwIW5ubowZM4ZGjRqxcuXKB+NwCSHo06cPH330EZqmcfu25Yu8hIQEvvrqqyyvVeRGkM9IcnIyw4cPZ/z48dSsWfOhbbt27eJ///sfXl5ehIWFMXjwYEJCivZI8ir7lKOyT0scKfu0JsrZkqO0RY7SFkuKg7YcP36chQsX0q9fvwfrTpw4wdmzZ9mxYwcdOnSgf//+NG/enO+++46IiAhcXFwIDg52nBHkM+Lh4cGCBQvYvHkzBw8efLA+ISGBjRs30q1bN1q1asXAgQPRNC3L4fmLAir7lKOyT0scKfu0JsrZkqO0RY7SFkuKg7Y0a9YMIQSLFi16sK5p06b079+fKVOmMHDgQP7880+++OILGjVqxEsvvfRQxUxGkXa2Lly4wOeff06bNm24dOkSVapUIS4ujkuXLuHk5MTLL7/8YF9N09i8eTN+fn588MEHBR1+gaCyTzkq+7SkOGSfeUE5W3KUtshR2mJJcdGW4cOHM2fOHBYuXIiPjw9ubm54e3tLTZvo6GjWrFnDkiVLHMvZMhqNfPvtt8yaNYtXXnmFwMBAGjZsyKxZsyhRogSDBg16qKIFpveu3bt3JyIiwkZR5x+VfcpR2aclxSH7zAvK2ZKjtEWO0hZLiou2zJ49m9jYWHbu3MnatWs5ffo0GzduZMmSJYSFhT2079KlS5k7d26W5yuSztakSZMoX748ZcqUwc/PDzC10SpVqhTNmjXL8rzLly8nODi4SL5SVNmnHJV9WlJcss/copwtOUpb5ChtsaS4acu1a9fw9vZ+0BMxPT2doKAgKlWqRFpaGt26dePcuXMYDAaGDRvmOM5Weno6N2/e5ObNm+zevZvly5cTGxtLWFgY586d49y5c1ke7+Pjw82bN7Pcx15R2acclX1aUlyyz9yinC05SlvkKG2xpLhpS+XKlR8a8sHV1ZVBgwbx/PPPP5iTuW7duhw9ejTL8xQ5Z+v//u//uHLlChcuXCA4OJjo6GimTZtGv379eOyxx3jzzTd57rnnKFmy5INjzp07R0BAAP7+/vz+++889dRTdOjQoTB/jlVQ2acclX1aUtyyz5yinC05SlvkKG2xRGnLP3zxxRcIIYiOjkYIwYwZMxzH2bp+/To7d+7k3XffBcDf358vv/ySxx9/HCEE3t7elChR4sH+er2eP/74g+3bt3Ps2DFSUlJo166drcLPFyr7lKOyT0uKW/aZU5SzJUdpixylLZYobfmHt956i7Zt2zJy5EimT5+e5b5FztkCk1NVt25d6bbhw4fTv3//B993795Nly5d2LhxIxcvXmTZsmUFEm9hoLJPOSr7tERln3KUsyVHaYscpS2WKG3JnHyNIG+PZFbR0uv1FrN4V6lShdDQUEaOHMnff/9dGOEVGCr7lKOyT0tU9ilHOVtylLbIUdpiidKWvFEkna3M0DSN+fPnc+TIEbp06UJAQAA3btxg7969fP/99+zcuZNnnnmmgCIueFT2KUdln5ao7FOOcrbkKG2Ro7TFEqUtmeNwzlZmCCEICgpi5syZnD17llmzZrFmzRqSk5NJSUkp0hUtUNlnZqjs0xKVfcpRzpYcpS1ylLZYorQlbziUsyVjwoQJ1KtXz2KQ06KIyj7lqOzTEpV9ylHOlhylLXKUtliitCVzio2zJePTTz91iIoWqOwzM1T2aYnKPuUoZ0uO0hY5SlssUdqSNxze2XIkVPYpR2WflqjsU45ytuQobZGjtMUSpS2ZU6ydLUdCZZ9yVPZpico+5ShnS47SFjlKWyxR2pI3lLNVhFDZpxyVfVqisk85ytmSo7RFjtIWS5S2ZI5ythwElX3KUdmnJSr7lKOcLTlKW+QobbFEaUveUM5WEUJln3JU9mmJyj7lKGdLjtIWOUpbLFHakjnK2XIQVPYpR2WflqjsU45ytuQobZGjtMUSpS15QzlbRQiVfcpR2aclKvuUo5wtOUpb5ChtsURpS+YoZ8tBUNmnHJV9WqKyTznK2ZKjtEWO0hZLlLbkDeVsFSFU9ilHZZ+WqOxTjnK25ChtkaO0xRKlLZmjnC0HQWWfclT2aYnKPuUoZ0uO0hY5SlssUdqSN5SzBcyYMYN3330XFxeXQrleXlHZpxyVfVqisk85ytmSo7RFjtIWS5S2ZI5ytrIgJiaGWbNmERcXZ+tQskVln3JU9mmJyj7lKGdLjtIWOUpbLFHakjcc1tlKTExk9OjRjBo1ijp16mS632effcbt27cZPnw4vr6+/PHHH/Tp0yfP1y1IVPYpR2WflqjsU45ytuQobZGjtMUSpS2ZYzNnSwjRVQhxXghxUQjxoXldAyHEASHEIiFEttcPCwtjyJAhvP/++wQHB3PgwAHu3buX5TH79+9nwIABpKWlZVnRSkpK4vr163h5efHNN9/wySefsH379tz+zEJDZZ9yVPZpico+5ShnS47SFjlKWyxR2pI3CszZEkI4A38DnYDrwBFgIDAa+BB4Ebigadr/MjleO3z4MEuXLqVbt264u7sTGRlJeHg4ERERpKWl4erqik6no3bt2jRp0oTatWszc+ZMUlNTiYiIYPbs2ZQoUSLTGBcvXszatWtZvXo1CxYs4M6dO+j1eiZMmGDt22EVVPYpR2WflqjsU45ytuQobZGjtMUSpS2ZYytn63HgoqZplzVNSwOWA70AZ0ADjIA0qPusWrUKX19fPDw8EEJQrlw5WrVqRe/evXnhhRfo06cP3bp1w9PTkx07djBhwgRq1aqFu7s7L730UpYVLYDKlStTvXp1du7cyalTp3Bzc7Pbihao7DMzVPZpico+5ShnS47SFjlKWyxR2iJn4cKFWW4vSGerH9BV07Q3zN8HAa2AH4F5wAXgVU3TDJkcr7355psIIRBC3F/3YHk/05Btd3JyYsaMGdnGGB8fz4QJEyhTpgzR0dGcPXuW1157jQEDBuT9hxcgKvuUo7JPS1T2KUc5W3KUtshR2mKJ0hY5Gzdu5LnnnrOf3oiaph3XNK2VpmmvZFbRus/gwYPp0aMHX3/9Nd27d+fLL7+ka9euTJkyhY4dO/LZZ5/Rvn17PvnkE9q1a8dHH33EU089xRdffMHWrVsxGAxZLj09Pbl8+TLXrl2jRYsWvPjiiyxYsIDTp0/n6PjCXs6bN4+ff/7Z5nHY23LMmDFs2rTJ5nHY03LZsmV8/fXXNo/D3pYzZ85k7dq1No/D3pbjx4/nzJkzNo/D3pZvvPEGKSkpNo/DnpaXL19mzJgxNo/D3pYNGjTIsu5TkM5Wa2CSpmldzN/HAWiaNiWHx9tnN0mFQqFQKBQKCZk5WwVZ2XLB1EC+I3ADUwP5lzRNCy2QCyoUCoVCoVDYIQU2ZLqmaXohxHBgC6ZG8QtVRUuhUCgUCkVxw24HNVUoFAqFQqFwBFQXC4VCoVAoFIoCRFW2FAqFQqFQKAqQQqlsCSEWCiFuCyHOZFg3XQhxTghxSgjxf0KIUhm2jTNP8XNeCNElw/oXhRB/CiHezbDuihDitBDihPlvVmH8pvwiuycZto0WQmhCiNLm70IIMct8T04JIZpn2HeU+Z68kGGdIcP9OHF/qqSiQGb3RQgxwlxeQoUQX2ZY7/Bl5VEy+f/UVAhx0Py7jgohHjevz1XZKSoUVDkRQrwjhJiZYfv3Qojtj5zfbstNQemKEOKbR/4vbRFC/JDh+1dCiPcK9tfljYJ6/jhiWRFC+AkhtgkhLpiXvub1xaKsFCiaphX4H9AOaA6cybCuM+Bi/jwNmGb+XB84CeiA6sAlwNm8bQ2mxvbLAW/zuitA6cL4HQV9T8zrK2PqVBB+/3cB3YHNmEbcfwI4ZF7vDfyKqaPD2gznSLD177NyWekAbAd05u9lilNZyeE92gp0y1Beduel7BSVv4IqJ0BL4HCGcx7E1JP6/v7LgBdt/ftzc1/M6/OlK0A/4DfzZyfgGHAgw/kPAE/Y+vfnoqzk+/njiGUF+BL40Pz5wwz3pViUlYL8KxRnS9O0vUDMI+u2apqmN389CFQyf+4FLNc0LVXTtDDgIqapf+Cf6X20DJ+LJLJ7YuYbYCym33ifXsBizcRBoJQQojwP3w+HIJP78hYwVdO0VPM+t83ri0VZeZRM7pEGlDR/9gEizJ8dsuwUYDk5AdQWQngIIXyAZPO6Rub9ngT2FcBPsgoFqCv7gdbmzw2AM0C8EMJXCKED6gF/Wu+XWI8CfP44YlnpBSwyf14EPJ9hvcOXlYLEXtpsvY6p1gxQEbiWYdt18zqA1cBR4KimafEZ9tmV4dXQqAKPtoAQQvQCbmiadvKRTdJ7Yr4HpzHdkxUZtnuIh18jFrnXRI9QG2grhDgkhNgjhHjMvL7YlhUJ7wLThRDXgBnAOPP63Jadoky+y4n5AXwceAxzBo/pYfykEKIiph7cGc9l91hDVzRNiwD0QogqmCoRBzDdm9aYHJ7TmmkO3KJInp4/jlhWgLKapt00f74FlDV/VmUlnxTYOFs5RQjxEaAHfsluX03TFvFPrTsjHTRNu2Pt2AoTIYQnMB6TvZ1jNNOI/I+Oyp+saVpTa8VmB7gAfpgE7THgNyFEjawOcOSykglvAaM0TQsRQgzANAczwCpGAAACjElEQVTps1kdkEnZKcpYq5zsx/SQ8MD0oLiA6f9mlHlbkcHKunL/vjwJfI3pAfwkcBc7dnCywgrPH4cpK4+iaZomcjCTS3EpK/nFps6WEGIw8BzwsqZp9/9Rb2BqX3CfSuZ1jk5NTG0ETgohrmD63X8KIcpRfO/Jfa4Dq80W9mHACJRG3ZeMvIop8wZYyT+vPorTPbJWOdmH6cHQGtMD9CymtjxPUvQeoNbUlfv3pRGmV0MHMd2jonhfrPX8caSyAhBpfj2IeXn/VXyxLivWwGaVLSFEV0xtCHpqmpaUYdM64EUhhE4IUR0IBA7bIsbCRNO005qmldE0rZqmadUwPTiaa5p2C9M9+be5R8gTwN0MVm9xYA2mxs8IIWoDbsAdimlZyYQI4Gnz52cwZdhQvMqOtcrJAUzuWICmabfND+IoTO1WilRWbmVd2Y+pchKjaZpB07QYoBSmh2iReoBa8fnjMGXFzDpMiRvm5doM64tlWbEauWlNn9c/TL0ybgLpmP6zD8HU8PAapgaFJ4B5Gfb/CFMvkPOYe1hlce4rmN4Z3z/P4sL4TQVxTyS/636vIQHMMd+T00DLbM5tyHA/TmBqNGzz35yPsuIGLMWUIf0JPFOcykoO71EbTL1+TmJqH9EiL2WnqPwVZDkx7x8KLM3wfRKQgLkHm73+FbCuOAP3gMkZ1v0MnLf1785DWbHK88fRygrgD+zAlKxtB/yKU1kpyD81XY9CoVAoFApFAWIvvREVCoVCoVAoHBJV2VIoFAqFQqEoQFRlS6FQKBQKhaIAUZUthUKhUCgUigJEVbYUCoVCoVAoChBV2VIoFAqFQqEoQFRlS6FQKBQKhaIAUZUthUKhUCgUigLk/wFRWVRMlORHYQAAAABJRU5ErkJggg==\n"
          },
          "metadata": {}
        }
      ],
      "source": [
        "# mo betta one\n",
        "if make_plots:\n",
        "    data = Dataset(traj_output)\n",
        "    traj = data.variables['trajectory'][:][:]\n",
        "    x = data.variables['lon'][:][:]\n",
        "    y = data.variables['lat'][:][:]\n",
        "    \n",
        "    fig, ax = plt.subplots(figsize=(10,10))\n",
        "\n",
        "    m = Basemap(projection='merc',llcrnrlon=120.0,llcrnrlat=0.0,urcrnrlon=280.0,urcrnrlat=60.0,\n",
        "                    resolution='c',lat_1=0.0,lat_2=60.0,lat_0=30.0,lon_0=180.0)\n",
        "    \n",
        "    m.drawcoastlines(color='black', linewidth=0.5)\n",
        "    m.fillcontinents(color='#c0c0c0', lake_color='#ffffff')\n",
        "    m.drawmapboundary()\n",
        "\n",
        "    m.drawparallels(np.arange(0,60,10.0),labels=[1,0,0,0])\n",
        "    m.drawmeridians(np.arange(120,280,20.0),labels=[0,0,0,1])\n",
        "    plt.title(model + ' ' + depth + 'm ' + lengthofrun + ' years')\n",
        "    \n",
        "# plot trajectory with a red start (r*) and the start/end\n",
        "# note the output (x,y) is dimensioned with trajectory number\n",
        "# first, then time.  The trajectory number alternates between\n",
        "# starting point and release time.  For example, x(0,0) is\n",
        "# the initial time, starting point 1, first release; x(1,0) is\n",
        "# initial time, starting point 2, first release; etc. then\n",
        "# second release and so on.\n",
        "    if repeat:\n",
        "        for i in range(0,int(ntraj*nstart),int(ntraj)):\n",
        "            for j in range(int(ntraj)):\n",
        "                traj_color = 'C' + str(j)\n",
        "                x2, y2 = m(x[i+j,:],y[i+j,:])\n",
        "                plt.scatter(x2,y2,facecolors='none',edgecolors=traj_color)\n",
        "        for i in range(ntraj):\n",
        "            x1, y1 = m(x[i,0],y[i,0])\n",
        "            plt.plot(x1,y1,'y*',markersize=16)\n",
        "    else:\n",
        "        for i in range(int(ntraj*nstart)):\n",
        "            x2, y2 = m(x[i,:],y[i,:])\n",
        "            plt.plot(x2,y2)\n",
        "        for i in range(ntraj):\n",
        "            x1, y1 = m(x[i,0],y[i,0])\n",
        "            plt.plot(x1,y1,'r.') # color + symbol... so r* = red star... bs=blue square\n",
        "    \n",
        "    plt.savefig('/content/drive/MyDrive/REU_2022_copy/plots_copy/Box/since_1980/' + experiment + '.png', format = 'png')"
      ]
    },
    {
      "cell_type": "code",
      "source": [
        ""
      ],
      "metadata": {
        "id": "T1WzDVc2fKRP"
      },
      "execution_count": null,
      "outputs": []
    }
  ],
  "metadata": {
    "colab": {
      "collapsed_sections": [],
      "name": "Parcels_trajectory_code",
      "provenance": []
    },
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.8.3"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}